{ "cells": [ { "cell_type": "code", "execution_count": 1, "id": "e54c8e89-35c1-4045-a4de-bc28fce771bd", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/Users/marcel/opt/anaconda3/envs/duplessis2021_JGR/lib/python3.9/site-packages/xarray/backends/cfgrib_.py:27: UserWarning: Failed to load cfgrib - most likely there is a problem accessing the ecCodes library. Try `import cfgrib` to get the full error message\n", " warnings.warn(\n" ] } ], "source": [ "import matplotlib.pyplot as plt\n", "import xarray as xr\n", "import cmocean.cm as cmo\n", "import numpy as np\n", "import my_functions as my\n", "import gsw\n", "from matplotlib.dates import date2num\n", "import matplotlib.dates as mdates\n", "import pandas as pd\n", "import glidertools as gt\n", "from sklearn.linear_model import LinearRegression\n", "import math\n", "\n", "import my_plot_params\n", "\n", "years = mdates.YearLocator() # every year\n", "months = mdates.MonthLocator(interval=1) # every month\n", "week = mdates.WeekdayLocator(byweekday=mdates.MO, interval=1)\n", "weeks = mdates.WeekdayLocator(byweekday=mdates.MO, interval=3)\n", "\n", "yearsFmt = mdates.DateFormatter(\"%d/%m\")\n", "\n", "mnthFmt = mdates.DateFormatter(\"%B\")\n", "\n", "lightblue = '#5499c7'\n", "blue = '#21618c'\n", "orange = '#f39c12'\n", "green = '#27ae60'\n", "red = '#cb4335'" ] }, { "cell_type": "code", "execution_count": 2, "id": "76a1bca0-cd5b-4766-b416-b154db2157da", "metadata": {}, "outputs": [], "source": [ "def rolling_mean(dat, window=4):\n", " dat_new = gt.cleaning.rolling_window(dat, func=np.mean, window=window)\n", " return dat_new" ] }, { "cell_type": "code", "execution_count": 3, "id": "43760fc1-7cc5-40df-9cc1-4274bfb2f891", "metadata": {}, "outputs": [], "source": [ "# note these are the 6H time series\n", "dat_saz = xr.open_dataset('../data/dat_saz_6h.nc')\n", "dat_pfz = xr.open_dataset('../data/dat_pfz_6h.nc')\n", "dat_miz = xr.open_dataset('../data/dat_miz_6h.nc')" ] }, { "cell_type": "code", "execution_count": 4, "id": "29e08cdd-5162-41ad-82ea-e0ed13b7802b", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAqAAAAFyCAYAAAAqImh8AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAB7JElEQVR4nO3dd3hUVfrA8e+Z9EZ6SAgJEHpHQZAmFlRs2HvHupafq6yu7qprXV1XdN21r2J37QXFiqKIgBTpHQIklPTey5zfH3cymTTSZubOTN7P88yTzM0tbyb3zbxz7j3nKK01QgghhBBCuIvF7ACEEEIIIUTPIgWoEEIIIYRwKylAhRBCCCGEW0kBKoQQQggh3EoKUCGEEEII4VZSgAohhBBCCLfyNzsAIUTrVs+c+DpwZTurPThh0coHXB9N16yeObE/sAc4f8KilR+ZHE4TaXPmnQDcBUwEQoC9wMfA4+nz55Z2Yj9XAa8B8enz5+Z1YjsN3Jk+f+6TnQjbJeRcE0K4mxSgQniuh4EXHZ6/Cey0LW+w360R+Yi0OfNOBb7AKBz/A1QARwD3AMelzZk3PX3+3HoTQ3Q3OdeEEG4lBagQHmrCopW7gd0Nz1fPnFgB5E5YtHKFeVH5jDuB79Lnz73WYdmPaXPmbQO+BE4GvjIlMhPIuSaEcDcpQIXwcqtnTrwY+AswGDgA/GvCopX/cfi5Bq4BTgNmAcUYLVsLgJeB4zBat26bsGjl17ZtfgK2AGXA9UAt8D7wpwmLVlYdJpZJwDzgSIyC5m+trDMG+AcwybZooW2/2bafvw6EAyuAPwLRGMXgNbbnNwN+wLvAHycsWmnt0AvVVAKtt+h9B/zV8Wdpc+ZNBB4ApgChGJd5n0qfP/eljh4sbc68IcC/gWlAlu13aL5Of+AJYIbtOD8Cf0qfP3en7ecPAKfb9nM/kAz8ClwOzLbFHYnRsnt9+vy5FR2Nr6PkXOvSuSaEaIV0QhLCi62eOfFKjDfHnzGKkDeAp1fPnHhns1WfxrikegawHHgWWIRRwJyPUSi8s3rmxFCHbS4Bjse4N/BB4Crg1cPE0h/4AagCzgPmY1zidlxnHMabfaBtv7cBxwA/r545Mcxh1ZOAc4DrMForzwFWA0fbtnsDuBW4sK142vE1cFLanHlfpM2Zd1HanHmJAOnz59amz5/79/T5czcApM2ZlwosxiiOzgfOBHYAL6bNmTemIwdKmzOvF/AT0Bu4FHjMFr/jOn2BlRiF3U3A1cAAYGnanHl9HFYdCvwZ4zW5FuP1+BmYY9vuMYy/222deTE6Qs61Lp9rQohWSAuoEF5q9cyJFuDvwDsTFq28xbb4O1sr1H2rZ058fsKileW25csmLFp5t227AxhvsssnLFr5d9uyKowiYQiwzraNPzBrwqKVebZ1NPDs6pkT/zph0cq9rYT0f0A1MHvCopUVwFerZ05UGK1UDe4DcoFTJixaWWPb7xpgI0YR1dCaFg6cO2HRykO2dS4HRgDjJyxaWQp8s3rmxAswWrb+18mXDozWwhiMAuN0ANvl948wWjcLbeuNxCiiLk2fP7fWtt5vQD5GMbOhA8e6CogHJqXPn5tp20chRoenBrdjdIQ6saEjU9qceT8B6cBc2wOM12VO+vy5v9nWOR24COifPn/uPuDLtDnzTqOxxc8p5Fzr1rkmhGiFtIAK4b2GAH2AhatnTvRveGC07kVg9O5usNLh+2zb19UOy/JtX6Mcli1qKAhsPrd9ndZGPFOBn20FQYOPm61zDPB5Q0EAMGHRyi0YhdwMh/UyGwoCh5i32woCx5gd4+2w9Plzq9Pnz50D9MNoOfwUo4XyXmBT2px5A2zrfZ0+f+5MwC9tzryxaXPmnYfRAgkQ1MHDTQU2NhSfNp8Djp2cjgEWO/ait33/A01fF03Tv1s2kGsrPht0+XU5DDnXnP+aCtGjSQEqhPeKtX19F+O+uYbHKtvyJId1WxtWqL17BA81e55r+xrTxvrRQPNhiLJaWSeblrKBXg7PuxJvp6XPn7s/ff7cF9Lnzz0H477QORitlQ8ApM2Z55c2Z96/gELgd4xWwIbXXXXwMC1eF1sP+9xm63TkdalopXe+01+XVsi5JoRwKrkEL4T3KrZ9vZmmrU4N9nRz/7HNnifYvua0sX6+wzpt7aMAo6WxuURga6ei66K0OfOOxmhhm91wKRsgff7cOuC1tDnzZgPDbYv/itEx5grgq/T5c8vT5swLxeik0lH5DvtriEFhFEgNDve65Ley3N3kXBNCOJW0gArhvbZhvBH3nbBo5eqGB8Yb8cMYPaK747hmHUXOAqzAkjbWX2zbJsph2anN1lkKnLl65sTAhgWrZ04cDozG6KTiDjswLhv/X/MfpM2Z5wekAZtsiyYDq9Pnz/0wff7chnscZ9m+drQFdDEwKm3OvMEOy06g6SX8pRjjj8Y5xBJnW89dr8vhyLkmhHAqaQEVwktNWLSybvXMiQ8AT62eORGM+wUHYPSE3olzWqUWrJ458SlgEMbl5+cnLFp5sI31/4XRWvj16pkTHwX6YruU7eBRYJltnacxCpdHMGYheoNuWD1z4hFAte0+vzalz59bkDZn3l+Bp2xF3usYQwP1AW6wxX2ObfVVwN1pc+bdgtF55SiMIZA0xlBJHfEm8CfgC9txQzBey1qHdZ7G6Kz0fdqceQ9jFLf3AjUYr2uXpc2ZNwIISp8/d21X9yHnWlMdPdeEEG2TFlAhvNiERSufBW7EGBbnK+Ah4EPgtAmLVupu7v5bjJbA9zE63vwTY2zEtmLJwejcUQ58ANyBUdA5rrMGY7idAFuczwC/AFObdfroik+B5zuyYvr8uU9jDBOkMcbV/BGjqMkEJqTPn9swKPvjGMXK3zAGqL8EY0ie7zFaRztyrCqM33krRrH7GEZxWeCwTiYwHTiIUbC+CuwDJqfPn9vdGYiex3htukXOtSY6fK4JIVqntO7u/w0hhK+xDQ5eNmHRytPNjkX4NjnXhOiZpAVUCCGEEEK4lRSgQgghhBDCreQSvBBCCCGEcCtpARVCCCGEEG4lBagQQgghhHCrThegSqnZSql2h7BQSk1RSi1WShUppQ4qpd5USrU2K4UQQgghhOhBOlWAKqWmAG/TzgwgSqnhGAMVlwIXYwzCPBX4VikV0LVQhRBCCCGEL+jQTEhKqSDgNowp18qBwMNvwS3AIeBcrXWtbR87MeYQPhFjEGMhhBBCCNEDdXQqzlOAe4A7MaZMm9vO+puBLQ3Fp81229cBnYpQCCGEEEL4lI4WoKuAAVrrIqXUA+2trLVubYqyM2xft3XwmEIIIYQQwgd1qADVWh/ozkGUUinAk8BqjDmXm//8euD6VjZ9WWv9cneOLYSQHBPClSS/hOi8Tg9Eb2sB/ZPWOryD66dgdEiKBKZorXd39FhxMTG6X9++nYrPbHV+HW1U9hxlZWWEh3foz+kx/OvrzA6h037fuDFPax1vdhyOJMfcQ3LMPTwtxyS/3MMb8wskx1x6pimlRgFfAwHAiZ0pPgH69e3Lr1994ZLYXKUgPMHsEDot69BBEpP6mB1Gp8SU5ZgdQqeFpPTfZ3YMzUmOuYfkmHt4Wo5JfrmHN+YXSI65bCB6pdQkYAlQD0zXWm9w1bGEEEIIIYT3cEkBqpTqj9HymY1x2X2nK44jnOOdN+abHYLwICEp/QlJ6c+S5cvNDsVnSI4J4TqSX97JKZfglVIDgXit9QrbomeAXsDNQKpSKtVh9X1a60POOK4QQgghhPA+zroH9D7gSkDZZjo6FfAD3m1l3TsxesQLIYQQQogeqNOX4LXWDzTvAa+1vkprrWzf12qtA7TWqo1Ht4vP4pISbrrrbtImTKJX2mD6HzmBG/90F4VFxd3ddY80edoxZocgPNCBQ1mcdcVVRA8ayoipx/DqO619nhQdITkmmnv+tdcZe+zxRA4cwqCJk7nj/gcoLSszOyyvJPnlnbxvvAXgujv+xBfffodSipioKLJz83jj/Q+wWCw8/8TjZofndaZI8opW3HL3X6itq8PPYmFPRga33P0XYqKjOfvUU8wOzetIjglHD/5zHo//+z8A9IqI4FB2Ni+89jq79+zh87feMDk67yP55Z1c1gveVWpqavD392fY4EFs/uVn9m9YyzOPPgzAqrXrzA3OS7347DNmhyA80LBBg9j3+2oObVrPScfOAOCxf8m50hWSY6JBQWERT7/0EgBP/O0+srds5NeFC/D39+fXlavYmZ5ucoTeR/LLO3ldC2hgYCDvvvg8VquVTdu28fKbb/HJwq8AKKsoNzk671ReVmp2CMID3XzN1URHRQJwxx9u4LuffmbLjp3U1dXh7+91/zpMJTkmGqxcu5bq6hqCg4K4ec7VAIwbNYpNv/xESp8+WCxe1y5kOskv7+SV7yKvv/c+D/7zSbJycukdH8eAfv0AsFqtJkcmhO9IiIuzf5/UuzcA9fX15OYXkNTb+warFsITFBYVARAVGdmk2PS2GZOE6C6v+6i1edt2brrrbgqLi/np80/Y+/tq/vrH28wOy6sl9E40OwThgTIOHLB/n52TC4BSisheEWaF5LUkx0SD6KgoAPILC6mpqbEv//jLhSz6eQlFxdKZtrMkv7yT1xWgW3fuRGuNRVlITkqiurqatz78CACrtXPz2gvD5Vdfa3YIwgP9+7+vcCg7h+rqap560bhnbdyokYSGhJgcmfeRHBMNJh5xBEFBgdTW1vKvl/4LwKat27jmj7dzxmVXsHXnLpMj9D6SX97J6wrQI0ePJiAggMqqKkZOm0Gf0eP44PMFAJTKfSBd8t3XC80OQXig8opKBk08mqRRY/nmx8UA3HPb/5kclXeSHBMNYqKjuOPGGwH42xP/pPeI0UyadSrV1TUcP30akyeMNzlC7yP55Z28rgBN69+P1/79LwYNGICfxUL/lBReevIJoiMjKS4p5fcNG80O0etsXL/W7BCEB/rglZc5deYJAKT168f8f/+LM04+yeSovJPkmHB0/5/u4MkH/8bgtDSqqqtJTkri5muu5t0Xnzc7NK8k+eWdvLIT0rmnn8a5p5/WZNkVF15gUjRC+JbKzL327z989b/mBSKED7t5ztX2XvBC9EReWYAKIYQQwvtVVVZQUV3TYnloUCDBIaEmRCTcRQpQwQ03yygCQriS5JgQrauoruHdJZtbLL/kmJEdLkAlv7yT190DKpwvK+uQ2SEI4dMkx4RwHckv7yQFqODzjz8wOwQhfJrkmBCuI/nlnaQAFUIIIYQQbiUFqBDCqQoKiygplTF5hRBCtE0KUMGJs041OwThoebcdju90gZzMCu7w9uMmXEcmQcOAvDEf55jzm23uyo8ryE5JoTrSH55JylABWPGHWl2CMIDFRYV880PP3LOaafyytvvdHi7/MJC+/d33Xoz85952hXheRXJMSFcR/LLO0kBKpj3+CNmhyA80Dsff8zUSRO54cormP/uu9TUGGP1PfLU01z9f3/knKvmEDd0BEccP5NFPy8BYMqppwMw/YwzWfDNtzzy1NNcfMMfTPsdPIXkmBCuI/nlnaQAFUK06rV33+PKCy9g8oTxxMfG8fHCr+w/+/jLhdx67TUc3LiOk487jjvufwCAZV99CcAvX3zO7FknmxG2EEIILyAFqBCiheWr11BUUsIpJxwPwLWXXcqLr79h//mkI4/guGlTCQwM5KKzz2TXnj1mhSqEEMILyUxIgrSBg80OQXiY+e+8S35hAQOPmgRAXV09+YWF/L5hIwBxsbH2dQP8A9BaH3Z/GQcOcOTxJ9qf/+fxv3Px2Wc5P3APJTkmhOtIfnknKUAFZ59/odkhCA9SXFLCJwu/4qv/vUtav1T78j/97UFeeO11+qX07fQ+U5OTydu+xZlhehXJMSFcR/LLO8kleMGnH75vdgjCg7z7yacM7N+PKUdNIDEhwf646qIL+fCLL8grKDzs9oGBgZSUlbkpWu8gOSaE60h+eScpQAXpu3eaHYLwIK+9+z/OP3N2i+XHT59GbHQMr/3vvcNuf8UF53HaxZfy9ocfuSpEryM5JoTrSH55J7kEL4RoYuV337S63GKxsHvVihbLRw4bSmXmXvvz/zz2d/7z2N9dFZ4QQggfIC2gQgghhBDCraQAFcy9+16zQxDCp0mOCeE6kl/eSQpQwYZ1v5sdghA+TXJMCNeR/PJOUoAKvv/mq/ZXEqKHyNi31+n7lBxrqaCwiJLSUrPDED5A8ss7SQEqhBA2mzas48yTj+vSttu2bCY5OrRL2y5e9B0XnHkqI9P6MnJAMpecO5v1a9d0aV/eYsyM48g8cLBL2/YdcwRLli93ckRCCHeSAlQIIWxKSkqora1z6zHfeWM+t998Pdf94RbWbd/Dmq27Oea4E7jgzFPZvtV3B+/PLzz8eLJCCN8mBajgzHMvMDsE4aGsVit//9czpE2YROLI0Zx/zXX2wuH3DRs58bwL6D1iNGOPPZ63PvjQvt3QyVN59tX5jJp+LLFDhnPrPX/h28WLGTltBr1HjObOBx6yrxuS0p8nn3+BfkdMIHn0OO5//AmsVisAlZVV3HH/A6RNmMSA8RO5++FHqampAWDe449w6/VzuOLCcxjcN55jjz6Sn39cBEBNTQ133HIDowamcOTwNK678hIKCvIB4xL7VRefx4SRgxmYFMPsk45j147t5OXmcPn5Z1FYkM/gvvEUFORTWVnJfX+ey/gRAzlyeBoP3XeP/fhWq5XHHrqfkWl9OXJ4Gp9/0vj7N9dWjlWUl/PQfffwz2ee58RZpxIQEEBwcDA33nIbV865np07tpObk83N117FqIEpTBg5mEfu/yvV1dUA/PGm63nk/r9y5snHMyg5jnNPO4m1a1Yx+6TjGNw3novPOYPSkhIAzjv9ZB65/69MmzCGISkJXHfFxRQWFthjeePVl5k6fjQj0/pyzWUXkpWTA8CS5cs56sRZ3PXgw/QZNZaBRx3NvBdetG/37/++wuBJU0gePY6Z515gn661srKK2/5yLyOnzSB2yHBGTT+WBd98C8CUU08HYPoZZ9qXvfzmW4yafizJo8dxwbXX248P8N6nnzN86nQSho/ir39/rN2pX0XPIu9h3qnTBahSarZSqt0bd5RSo5RSPyilypRSGUqpPyulVNfCFK6UmJhkdgjCQ736zru8/eHHfPP+u+z7fTVhoaHccd/fyM3P59SLL+Hs005l//rf+e/T87j74Uf5dvFi+7afLPyKX774nGVffcFr/3ufp154iV8XLuD7j97nxTfeZMv2HfZ1v170I2t++I5fvvicDxYs4NV33gXgnkceZceu3az6/htWfvc1v2/YwL/n/cO+3Reffcx1N93K5vQDHH/iydz757kAfPz+u+zYvo2VG7bx65pNVJaX8+qLzwHwp/+7iUGDh7Ji/VY27sokNi6OZ+Y9QVx8Am99+BnRMbHs3J9LTEwsD993D7t27uD7pSv5/pff2LD2d/vx33z1ZRYu+Ixvf1rGj8vXsGbVyjZfx7ZybNVvy6mrq+O4mSe1+NlfHniY0888m2suuwilFCvWbeGL739i+a9LmPfYI/b13n/3Lf75zHOs376X3Nwc5lx6AU89+yKrNu1gf2YGH7//rn3dj95/h5dff4fft+ymurqav8z9o+11/IRnn36S+W+/z5otu0jtP4DLb7rFvt2mbduIiYoiY90annroAe5//An2HzpE+t59PPTkUyz6+AP2b1jLsVMnc9eDDwPwr5deZtuuXSz76gtytm7iigvO5477HwBg2VdfAvDLF58ze9bJfPzlQv753At88MrL7F61ggGpqfbjb9y6lZvu+jMvPfkE+9f/jlKKgqKiNl9r0fPIe5h36lQBqpSaArwNHLaQVEolAIsADVwAvAw8CsztWpjClV567hmzQxAe6oPPF3DTnKsYMnAgQUFBPPng3/jzrbew8PtF9E3qw01XX0VAQAATjziCOZdezNsffmzfds7FFxEdFcnQQYPsU3lGRUYyZsQIEhMSyDhwwL7uw/f8mbiYGNL69+PmOVfzwecL0Frz5gcf8shf/kxsdDTxsbHcd8ftvPvma/btxh81iekzjiMwMJBzzr+QPbt3ARDRK5I9u3fzwf/epiA/nzc/+JQ7/3I/AE8/9zJz77mXuro69mdmEB0TQ9ahlvciaq15/923+OsDDxMTE0tsXDxz77nXfvwFn37M1dfdSN/UVKKiornznvvafB3byrHCggKioqLw9299TpC9e9JZs+o3Hnr8n4RHRJDUJ5k7/3I/H/zvLfs6M0+exZBhwwkLD2fsEeM54aRZDBoylKioaI4YfxT7MzPt61593R8YMWo04RER3HXv3/hm4QKqq6t57+3Xue6mWxk6fATBwcHcc/9DrFq3jp3p6QD4+fkx96Yb8ff358xTZhEeFsbefRmEhoZQU1vL/Hf+x4YtW/jLH29j0ccfAHDDlVfw7ksvEB4Wxv6DBwkPD+NgVlarv+cb773PrdfOYcTQIQQHB/Pw3XfZj//pwq858dgZHDN5MoGBgdw/9w7CQrt2r63wTfIe5p06NBOSUioIuA14GCgHAtvZ5GbbvmdrrSuAr2z7uEcp9YzWurYbMQsh3CQnL4/kpMbWhbiYGOJiYlj4/fek9k1usm5qcjK/rlxlfx4dFWX/3s/PQmSvXvbnFouyX2YHGNi/n/375KQksnPzyM3Pp7KqipMvuJiGaydaQ01tLVVVVQDExMXZt/P3D7Bfmj39zLPJz8vlg3ff4v67/8SwESP5x9P/4YjxR7F753auvuSvZB06xJBhw1FKoR1iaZCfl0tVZSXnnT6Lhos3Wmtqa2uoqqoiJyebxKQ+9vX7pqa2/4I2E9+7N0WFhdTW1hIQENDkZ0VFhWQfOkRoWBgxsY2/Z9+UVHJzcqitNf6NRkXH2H/m5+dHr8go+3OLxYJVN/5u/dMG2r/v0yeZmpoaigoLOLB/P088+iBP/6NxBiuFImP/AQIC/Inq1atJfAH+/li1lcSEBD5/83WefvEl/vPKq0RHRfG3P93BFRdeQElpKbf99V5WrV3HgH79GJCa0ual88yDB3nwn/P4+78aC4mG42fn5tInMdG+PDAwkMSEhHZfWyGEZ+voVJynAPcAdwKxtN+SORP4wVZ8NvgMuBc4CljWuTCFEGbok5jYpNVqb0Ym73z8MWn9+7Pv08+brLs3M5MEh4KwM3fcHMrOpnd8PAAZ+w/Qt08SsdHRBAYGsuLrhQzoZxR35RUVbC+3EhwcfNj9pe/exdRjjuXKa66noCCffz3xGLf94ToWLV3JtVdczFPPvsTpZ54NwNNP/J1fl/zcYh/RMbEEBgby7ZLl9Os/ADDu2czJySY4OJjExCT2Z2bY12+tFbU944+aREBAIIu//5aTbPdFNvjTrX+grKyMivJyCvLz7EVoxr69REXH2AvCzrzO2VmH7N/vz8wgOCSE6JhYEnoncuMtt3HRZVfaf56zdhlp/VJZsabt3vi5+fmEhYWy4O03qaqq4tOvvuaaP97BzBkzuOXuvzBsyCA+fu1V/P39WbriNz7+cmGr+0lMSOCP11/PlRc13su3becu0vqlsnz1atZt3GxfXldXR25+fod/ZyGEZ+roJfhVwACt9b8xLqu3Zwiwq9mydIefCQ8yeuwRZocgPNRFZ53FC6+9QfrefVRVVfHQvKfYtWcvs447jpy8PJ5/7XVqa2tZuXYtr737HhedfVaXjvPIvH9RWlbGzvR0nn/tdS4771z8/Py46Kwzuffxf1BUXEx5RQW33P0Xbr/p+nb3991XX3LztVeSm5NNVFQ0oWFhREfHUFtTQ1VlJaG2S7hrVq3krddesbcmBgYGUV1dRU1NDX5+fpx9/oX8/cH7KC4uoqK8nD/ffqv9+OdeeDGvvPAsu3buoLSkpMl9mc21lWPBwcHcff+D3HX7LSz69mvq6uooKy3l6Sf+zi8/L+aBR//BtBnHcf/dd1JeVsahgwd48rGHOef8Czv7EgPw+n9fZN/ePZQUF/PEIw9y1rnnExgYyPkXX8pLz/2bPem7sVqtzH/5BWaceTblFZWH3V/G/gOcfsnlrN24ieDgYGKjowkOCiIsNISSslJCgoPx8/MzWjifnAfg8FoHUlJWBsBl553LM//9L7v37MVqtfL8a6/bj3/+7DP44Zdf+GrRD9TW1vL3f/1bxg8VTch7mHfqUAuo1vpA+2s10Qto/h+i1OFnwoOcdMppZocgPNQVF55PTl4up1x8CaWlZZxwzHT+89ijREVG8vlbb3DnAw/xwBNPEhcTw8P3/JmzTpnVpeOkpvTliONPpL6+jluuvYZLzj0HgCcf/Bv3/v1xjjzhJCqrKply1FG8MP/Ndvd3zY03s3dPOidMnUhVVSVjxh7BU8+9RFh4OI8/9W/uvO0mysvK6TdgAJdddQ2vv/IydXV1jBg5iqHDhjNqYF++/Xk5Dz32JI8+cC/HTx5PZUUlRx092X78iy67kuzsLM459US01lwx51oW//B9q/EcLseuuvYGIiOjePqJv/N/N1yDslg4cvwEPvriG4aNGMmzL8/n/rvv5OhxIwA454KLuOf+h9rc3+GMP2oSV198PocOHmDW6Wfw0GNPAnDehZdQVFjIZeefRV5uDoMGD+GT1+cTHRV5+P2NHcPDd/+Zi2+4kdy8fFL7JvP2C88R2asXT/ztPm7581948fU3iYuN4drLLmXtxk1s27WL0cOHc8UF53HaxZfy778/wqXnnUtBURFnXnEVOXl5DBk40H786KhI3nzuP/z5oUe4Ijub82efwcD+/bv0+wvfJO9h3kl1djgLpdQDwJ+01uGHWacGuE9r/Q+HZf5ALXCr1vrZZutfD7Ro1khJTh6/Y8WvnYrPbAXh3ndv0luvvcLlV19rdhidElOW0/5KHiYkpf8arfUEM47tyTkWktKf1d9/y8hhQzu0vuRY15x3+smcNvssrr7+Dx1aX3Ks4zw5vzrL3flVUFTEu0s2t1h+yTEjiXG4j/xwPCG/uqKn51hH7wHtrGIgotmyCIefNaG1fhmjp3wT48eMkcHe3CAnu/WeqcJ3SI6ZS3LMt0l+mUvyyzu5aiD6nUBas2UNz7e76JhCCCGEEMILuKoF9AfgBqVUmNa63LbsLCAfWOeiY4ouCgtv3lgthPtUZu41OwSX84Qc++jLb80OQQiX8IT8Ep3nlBZQpdRApdTRDouexxgr9Cul1OlKqb9iDOP0uNa6xhnHFM5z4y23mR2CED5NckwI15H88k7OugR/H7C84YnW+hDGWKD+wEcYN2f/VWv9pJOOJ5xo2dIlZocghE+THBPCdSS/vFOnC1Ct9QPNe8Brra/SWqtmy1ZrradqrYO11v0ce8QLz7JcklcIl5IcE8J1JL+8k6s6IQkhhBBCCNEqKUCFEEIIIYRbSQEquPTKOWaHIIRPkxwTwnUkv7yTFKBCCCGEEMKtpAAVvPPGfLNDEMKnSY4J4TqSX95JClAhhBBCCOFWUoAKIYQQQgi3kgJUMHnaMWaHIIRPkxwTwnUkv7yTFKCCKZK8QriU5JgQriP55Z2kABW8+OwzZocghE+THBPCdSS/vJMUoILyslKzQxDCp0mOCeE6kl/eSQpQIYQQQgjhVlKAChJ6J5odghA+TXJMCNeR/PJO/mYHIMx3+dXXmh2CED5Nckz0dFWVFVRU17RYXl9X3+19S355J2kBFXz39UKzQxDCp0mOiZ6uorqGd5dsbvGot1q7vW/JL+8kBahg4/q1ZocghE+THBPCdSS/vJMUoEIIIYQQwq2kABVCCCGEEG4lBajghptvMzsEIXya5JgQriP55Z2kABVkZR0yOwQhfJrkmBCuI/nlnaQAFXz+8Qdmh2CagnoLGZW6xaOgXlJDOE9PzjEhnKmqsoKCoqImj6yD+80OS3SBjAMqerSymnre+3lDi+UXzRhDTIgyISIhhBBtaRjOydGEhBCTohHdIc08QgghhBDCraQAFZw461SzQxDCp0mOCeE6vRNlKk5vJAWoYMy4I80OQQifJjkmhOtERkWbHYLoAilABfMef8TsEITwaZJjQrjOjm1bzQ5BdIEUoEIIIYQQwq2kABVCCCGEEG4lBaggbeBgs0MQwqdJjgnhOuFh4WaHILpAxgEVnH3+hWaH4PUK6i2U1dQ3WRYe6EeMn9WkiIQnkRwTwnX6pKSYHYLoAmkBFXz64ftmh+D1Gga0d3w0L0hFzyU5JoTrHMzMNDsE0QVSgArSd+80OwQhfJrkmBCuU1ZeZnYIogs6XIAqpa5TSu1USlUqpZYrpSa3s/4UpdRSpVSpUipdKfU3pVRA90MWQgghhBDerEMFqFLqCuBF4G3gXKAI+FYpNaCN9QcC3wFltvWfBv4MPNb9kIUQQgghhDdrtxOSUkoBDwEva60ftC37HtgO3A78XyubnQf4AedqrcuB75RSScAtSqk7tdbaWb+A6L65d99rdghC+DTJMSFcZ8iw4WaHILqgIy2gg4B+wIKGBVrrWmAhMKuNbYKAWqDSYVk+EG77mfAgG9b9bnYIQvg0yTEhXKe4qNDsEEQXdKQAHWL7uqvZ8nRgoFLKr5Vt3gHqgceUUjFKqaOAPwKfaq2ruhqscI3vv/nK7BCE8GmSY0K4TnZWltkhiC7oSAHay/a1tNnyUtv2Yc030FrvBv5ke+QDK4Ec4OouRyqEEEIIIXxCRwaiV7avze/bbFjeYqRtpdS1wH+Bl4H3gT4Y95EuVErN1FpXN1v/euD65vtJSU7uQHhCiPZIjgnhOpJfQnReRwrQYtvXCCDbYXk4RvFZ3so2dwNfaa1vaFiglFoNbAUuBeY7rqy1fhmjWG1i/Jgx0lnJDc489wKzQxAuJjlmLskx3yb5Za4+yTITkjfqyCX4hhGU05otTwO2t9GjPQVY4bhAa70N43L8iM4GKVwrMTHJ7BCE8GmSY0K4TnBwsNkhiC7oaAGaCZzVsMA2oPxpwA9tbLMDmOq4QCk1CIgF9nQlUOE6Lz33jNkhCOHTJMeEcB2Zacw7tXsJXmutlVKPA88qpQqBX4FbgDiMAeYbBp6P11o3tHo+BHyglHoF+B+QCDwA7AXedPLvIIQQQgghvEiHZkLSWj8P3AlcDnwERAEna63TbavcByx3WP9DjBmQjgS+wpgBaQkwSWvdvDe9EEIIIYToQTrSCQkArfU8YF4bP7sKuKrZsk+AT7oRm3CT0WOPMDsEIXya5JgQrhMZFWV2CKILOtQCKnzbSaecZnYIQvg0yTEhXKe3dPLzSlKACt567RWzQxDCp0mOCeE6+/ZI32Zv1OFL8MJ35WTLNGZCuJLkmBDtq6iuZeOeLLIKSln421bGDUzmyhPHMywl4bDbVVfLDN/eSApQIYQQQpgqv6SCXzbuoaq2zr7sl017WL51H7edNY2zp44yMTrhCnIJXhAWHmF2CEL4NMkxIdpWWV3bovhsUFdvZd7HS/j+97bH+vT3D3BleMJFpAAV3HjLbWaHIIRPkxwTonVaa1Zsy7AXn4H+fkwb1Z95157M0L7x9vUee+9H0g/lt7qPtEGD3BKrcC4pQAXLli4xOwQhfJrkmBCt+313FtmFZfbnU0b0o29cJIP7xPL0jWfQLyEKgJq6eh57fzH1VmuLfeTn5borXOFEUoAKlsuboxAuJTkmREtaa979eaP9+aA+sSTGNN6u0is0mEeumkWAn1GqbM3I4Zs1u1rsJz8vz/XBCqeTAlQIIYQQbncwv4R9ucUA+FssjOrfu8U6AxJjuPKkCfbn7/28kdq6erfFKFxHClAhhBBCuN2ug433dA5KjiU4sPXORBcfO47e0eEAFFdUsy1TLrn7AilABZdeOcfsEITwaZJjQjRVWlnNoYJS+/NBfWLbXDcowJ9rZ020P995IK9JK2hqvwGuCVK4lBSgQgghhHCrPYcK7N/3iYkgPCTosOufNH4IfWJ7AUaHpPSsgsOufzhVlRUUFBW1eFRVVnR5n6LzpAAVvPPGfLNDEMKnSY4J0UhrTUZOkf35gKSYdrfxs1i46Nhx9ufbM3OxWjUAGfs6NxVnRXUN7y7Z3OJRUV3Tqf2I7pECVAghhBBuU1haSVmVUeyFBPrTJ6ZXh7Y79aih9Ao1WkorqmvJyC1yVYjCDaQAFUIIIYTbOBaOE4ck4+fXsVIkODCA048aYn++LSMHrbWzwxNuIgWoYPK0Y8wOQQifJjkmhEFrzf68YvvzKcNSOrX9qUcNwc9ilC5F5VXkFJUTGxfn1BiFe0gBKpgib45CuJTkmBCG0spqyiqNy+/+Fguj+yd0avteoUEMSIy2P991MI/YuPjDbCE8lRSggheffcbsEITwaZJjQhgO5pfYv+8dE06gv1+n9+E4ZNP+vGLWbdzilNiEe/mbHYAwX3lZafsrCSG6THJMCINjAdowrFJrlNYUFBW1WF5fV09UeAhxvULJK6lAa1i2M5vjp7siWuFKUoAK0QkF9RbKalpOAyeDdwghxOHV1tWTW1xuf3643u/VtXV8tGJHi+XnHW10QhrUJ468kgwAVqTnUW+12u8NFd5BClBBQu9Es0PwGmU19bz384YWy2dPG2NCNMJbSI4JAbnF5TR0Wo8KCyYkqPWpNzsiJT6S33f5UVNXT1FFLcu3ZjBtZH/nBCrcQj4uCC6/+lqzQxDCp0mOCQHZhY23ovSOjujWvvz8LKQ5DGD/2a+burU/4X5SgAq++3qh2SEI4dMkx4SArMIy+/eJ0eHd3t9Ah85Iv23PaHJ/qfB8UoAKNq5fa3YIQvg0yTHR0xWWVVJcXgWARSnio8K6vc+IkCCSYoyWVK3h8+Wbu71P4T5SgAohhBDCpbZk5Nq/j+0Vir9f54dfao3jkEwLf9tKTV3LTqLCM0knJCGEEEK41JaMHPv38ZHdb/1skBTbi6jQAIoqaikqr+Kn9bs5abzRU76qsoKK6pZjlNRLkeoRpAAV3HDzbWaHIIRPkxwTPd2WzMYWUGcWoBalOHXCUN5dYnRC+vTXTfYCtKK6hneXtLws3zCUkzCXXIIXZGUdMjsEIXya5JjoycqratiTVQSAAuKcWIACTBvWB38/o5zZuDeLnQfynLp/4RpSgAo+//gDs0MQwqdJjomebNPeLKy2AUCjwkMI6ML0m4dTUZTLjDFp9uefLZMhmbyBFKBCCCGEcJn16Y1XAJzd+tngnKmj7N9/t2YHZZXVLjmOcB4pQIUQQgjhMhscClBn3v/paMyAJNISjYHpK2vq+HZ1y2k8hWeRAlRw4qxTzQ5BCJ8mOSZ6qpq6erZkZNufu6IA7Z2YiFKKsx1aQT9dtgndMO+n8EgdLkCVUtcppXYqpSqVUsuVUpPbWT9eKfWmUqpAKVWklFqglEo73DbCHGPGHWl2CEL4NMkx0VNty8yxj80ZHhLYrfnf2xIZFQ3AyeOH2Pe/N7uQjXtzDreZMFmHClCl1BXAi8DbwLlAEfCtUmpAG+sHAN8DE4HrgKuAgcDXSqnAbkctnGre44+YHYIQPk1yTPRU7rj8vmPbVgBCgwOZNWGoffmC37a55HjCOdotQJVSCngIeFlr/aDW+itgNpAH3N7GZlcAQ4CTtNYfa60/Ay4FIoDRzghcCCGEEJ5t877Gy++u6oDk6LxpjSXGqh0HKKmQzkieqiMtoIOAfsCChgVa61pgITCrjW3OBr7RWmc4bLNOa91Ha72mG/EKIYQQwgtordniWID2cn0B2q93NFNG9DOOD+zYn3v4DYRpOlKANkwZsKvZ8nRgoFKqtQG9xgDblFJ/U0plKaWqlVILlVKp3QlWuEbawMFmh+BxrBYLGZW6xaPlpG5CtE9yrKWC+tZzrKBe+sb6iuzCMvJLKwAICfQnIjTIJccJDwtv8vzCGWPt3+/JKqC6ts4lxxXd05GpOHvZvpY2W16KUcCGASXNfhYPXA3sBa6xrfMPYKFS6gittZwNHuTs8y80OwSPU1VrZcHSDS2Wz542xoRohLeTHGuprKae935umWMXzRhDTIgyISLhbJsder8P7hOLRbnm79onJaXJ8yMHJTOoTyy7DuZTb9XsPpjPiH69XXJs0XUd+ajZcMY0H8+gYbm1lW0CgEDgFK31Qq31B8D5wCjgnBYHUOp6pdTq5o/cgoKO/RaiWz798H2zQxAuJjlmLskx3yb51TrHy+9D+8a57DgHMzObPFdKNWkF3XEgj3pra6WKMFNHCtBi29eIZsvDMYrP8la2KQN+01oXNSzQWq/G6D3fohOS1vplrfWE5o/4mJgOhCe6K333TrNDEC4mOWYuyTHfJvnVOscCdEhyrMuOU1Ze1mLZzCMGExMeAkBVTR0ZOUUuO77omo4UoA3/OZuP4ZkGbNetj/S6C6MFtDl/WrakCiGEEMKH1NXXs92hA5ArC9DWBPj7cdrEIfbn2zNzZWB6D9PRAjQTOKthgW2cz9OAH9rY5jtgqlKqj8M2MzBaTZd1NVghhBBCeL5dB/PtA9AnRkcQbWuNdKdZ4wfhZzHKnKLyKrILW7aUCvO0W4DaWjgfB25USj2qlDoV+ByIA54GUEoNVEod7bDZ0xiX7r9WSp2llLoEeBej+PzOyb+D6Ka5d99rdghC+DTJMdHTbMlonIVopIs7AA0ZNrzV5REhQaQlRbcakzBfh8a70Fo/D9wJXA58BEQBJ2ut022r3Acsd1g/F5gK7AHeAp7FmBnpNK213AnsYTas+93sEITwaZJjoqfZvC/L/r2re6AXFxW2+bOhfRPsPaZzisooKKlwaSyi4zo84JrWep7WOlVrHaq1nqK1diw4r9Jaq2br79Zan6W1jtBax9jWKXJi7MJJvv/mK7NDEMKnSY6JnmbLvsbWxhGpCS49VnZWVps/Cw8JJDUhyv58a6a0gnqKjowD6jEK6i2U1dS3WO4f4E9dKwPNhgf6EePXtMG1rX20tq4QQojuaet/bluTOhiTQDT9Xyz/n71LSUUVmblFAPhZLAxJjqe8wrz7L4elJrDP1gs+M7eYQwXNhzUXZvCqArStgYtnTxvT6qDhrQ1oLIMfCyGE+xzu/3ZrWpsEQv4/e5etDvdaDk6OJSjQn3ITr3xHh4eQFBNhLzwXrNxBbFTzkSWFu8mcZ4Izz73A7BCE8GmSY6In2eww/ueIVNfPQNQnOaXddYalNN4GsHjjXiqra10ZkugAKUAFiYlJZocghE+THBM9ieMA9O6YAjM4OLjddRKiwoiJCAWgrt7KDocxSoU5pAAVvPTcM2aHIIRPkxwTPYXWuklHH3e0gHZkpjGlVJPOUI7jlApzSAEqhBBCCKc4kFdCcXkVYIzDmRIfaXJEjZLjehEREgRAbb2V3QfzTY6oZ5MCVAghhBBO0Xz8T6U8p/OYUophqfH259v351JvldEVzCIFqGD02CPMDkEInyY5JnoKx9mGXD3+Z4PIqKgOr9u/dzQx4cY9o1U1dezNansQe+FaUoAKTjrlNLNDEMKnSY6JnsLdHZAAeneik5+fxcJpRw2xP9+WmYtVa1eEJdohBajgrddeMTsEIXya5JjoCapr69h5MM/+3B0dkAD27dnTqfVPHJdGgJ9R/pRWVnMgr9gVYYl2eNVA9MI1crLbnsZMCNF9zsixqsoKKqpbzh8UGhRIcEhot/cvRHftPJBHXb1xT2Xf+Egiw9ofHskZqqurOrV+aFAAg5Lj7APmb83IoW+c53SW6imkABVCCC9QUV3Du0s2t1h+yTEjpQAVHsHx8vtIN7V+dtWQ5Di22y6/F5RWklNUbnZIPY5cgheEhcuUZEK4kuSY6Ak2OfaAd2MB6u8f0OltQoICGJAYY3/uOH2ocA8pQAU33nKb2SEI4dMkx0RP4DgF58j+7itA0wYN6tJ2w1LiaRgkKquwlHTpEe9WUoAKli1dYnYIQvg0yTHh6/KKy8kuLAMgKMCfQX1i3Xbs/LyuTasZERpEX4eB8j/+dYuzQhIdIAWoYLm8OQrhUpJjwtc5DkA/PCUBfz8/tx07Py+v/ZXaMNxhrNJft2RwIF96xLuLFKBCCCGE6JZNe90//qczxESE0jsqHACr1rz303qTI+o5fLoAtVosZFTqJo+Wg5gIIYQQojsc7/8c5cb7P53BsRV04W9bKSytMDGansOnh2GqqrWyYOmGJstmTxtjUjSe69Ir55gdghA+TXJM+LLaunq2ZTb2Ih/ZL9Gtx0/tN6Bb2/eODic6PITCskpq6ur5aOlGrjtlkpOiE23x6QJUCNE5BfUWymrqmywLD/Qjxs9qUkTOJwO6u0Zr5w4gV516gN2H8qmpM/72STERxPbyrjxSSjE8NYFlW/YB8MnSTVx2/JGEBHV+eCfRcVKACt55Yz5z777X7DCEByirqee9n5teNbhoxhhiQlQbW3gfMwZ07wk51tq5A3LVqSfYtLexA5K7Wz8BMvbtoU9S947bNz6SxOhwsgrLKK2s5oe1Ozn96BFOilC0xqfvARVCCCGEazl2QBrpRR2QHFmU4pQJg+3PP1m2Ga21iRH5PilAhRBCCNFl3twBydEJY9MI9DeGj9qxP1dmR3IxKUAFk6cdY3YIQvg0yTHhzaoqKygoKmrxqKqsoKC0gkMFJQAE+vsxqE+c2+OLjXPOMXuFBjHziMZW0E+XbXLKfkXr5B5QwRR5cxTCpSTHhDc73H3Tm/c3Tl85NCWeAH/3DUDfIDYu3mn7OmvKSL5atQ2AH9bu4pbZU4kMC3ba/kUjaQEVvPjsM2aHIIRPkxwTvsqxA9IoEzogAaTv2uW0fQ1PTWBoX6OgramrZ+HKrU7bt2hKClBBeVmp2SE4VUF9ywkIMio1BfVyugtz+FqOCdHAE+7/rKurddq+lFKcPXWU/fnnyzZjtUpnJFeQS/DC57Q1HIyvDSckhBBmqrdaTR2A3lVmHjGIZxcso6yymgP5JazakcmkYalmh+VzpAAVJPT2jX8aXaW1pqiyhoNFFZRV1VJbb6WgFgoLS+gTFUpIoKSJ6J6enmPCN6VnFVJVUwcYswnFRYaZEkdQkHPv0QwODOC0icN4/2djXvjPlm2WAtQF5J1VcPnV15odginqrZrFmzP4fH0GRRVN52vZlWtcMlVAv9hwxqbEEB0aZEKUwhf01BwTvm3zvsbWz7EDkkyLo9+A7k3F2ZozJ4+0F6DLtuwlv6TC62Z48nRyU5zgu68Xmh2C2+3KKeGS//7IPxasalF8OtLA3vwyFqzPYH1mAfVW35mSUrhPT8wx4fu2ZOTavx+T1se0OLKzDjl9n6kJUYxNM4rqeqvm06XrmgxBVVZS3ObQVKJjpAVUsHH9Wk465TSzw3Cbz9ft428L1lBT11hM+lsUydFhxIcHE+hvISUxlsUb9pBTWgWA1rA2M59HP17GyIRwAvzks5vouJ6WY8L3aa3Z7FCAjhtoXgtocVGRS/Z72sThrE83ituPf91KZa0VpYx+BOcdPYSPVuxosY0rp/T1NVKAih5Da81T32/ilV+225cF+FkYnhTJyD7RBDmMXzd76ghCdB15ZVWs3JNrL0TX7ckmM6eIk0YmN1lfCCF6kpKKakorqwGICgumX0K0yRE533FjB/L0J0uorKmjtLKavJIK4k26z9UXdbgZRyl1nVJqp1KqUim1XCk1uRPbPqCUknEMhGm01rzw/fomxeeghF785+oTODI1rs1iMi48mFmj+jKmb+M/1/zyan7cdpC6erkcL4TomXKLyuzfj0lLsrcM+pKQoACmj+pnf55+KN/EaHxPhwpQpdQVwIvA28C5QBHwrVKq3Tt/lVKjgHu6EaNwsRtuvs3sEFxKa83qvXksWLPbvuzYIUn877rj6B/fq93tLUpxZGocR6c1zraRXVLF0l3ZaC2fq0T7fD3HRM+TW1xu/36sifd/AqQNHNz+Sl104riB9u8zc4qprat32bF6mnYLUGV8rHkIeFlr/aDW+itgNpAH3N7Otn7Aq0Du4dYT5spywQ3cnmTTwUI2HyqyPz9lVF+euXgyYUEBndrPsMQorjh2tP353vwytjjsV4i2+HqOiZ5Fa01OkwLUvPs/Aaqqqly27yHJsfSyjYBSZ7WSkVPksmP1NB1pAR0E9AMWNCzQWtcCC4FZ7Wx7O9AL+E9XAxSu9/nHH5gdgsvsyy9jzb7GyyYzh/fhH+dO7HInorMmDmFo70j789V788gtdd0/P+EbfDnHRM9TXlVLZbUx+1BIUACD+sSZGs/BA5ku27dSirSkWPvz9KwClx2rp+lIJ6Qhtq/NJ1tNBwYqpfy01i3apJVSg4AHMIrUCd0JUoiuKCyv5pedjfMUj06J45/nTcK/mz3YJw6II7+8iryyajTwy64sLquVyzLeTmlNQSu9aUODArvdq3XS5Cmt7jvQoqhpNs2fM47na6wWCxmVLe+5Dg/0I8ZP7sV2t9zixvs/R/dP7Pb/1O5KSIhvNb/qnXS5vH/vaDakH8KqNfklFRSXS6ODM3SkAG24Sa75ZMalGC2oYUCJ4w9sl+1fAd7SWi9VSkkBKtyqtt7KTzsOUWd7c48IDuC+c44mKKD7Pdf9LBZmDEni83X7qLNqSiprefeXzcQGdnvXwkTVtXUuG1bFPyiEd5dsbrG8taFcZBiXlqpqrSxYKtPreoqcZh2QzGbVqs38cobgQH/6xPZif14xAOmHpBXUGTpSgDZkd/PeFg3LW/v4eQPGpfvZHQlCKXU9cH3z5SnJyR3ZXHTTibNONTsEp1uRnkNxpXGJyN+iOGFYkv0+HmeICA5g4oB4lu02ZgJZuGYnp49JJSbMM2dLkhwzV+/ERMg5aHYYwkV6Un5prckubCxAJwzua2I0hqBg507F2Zq0pBh7Abo3u5BaGQWl2zrSbl5s+xrRbHk4RvFZ7rhQKZUCPAHcBlQopfwbjqOU8ldKtTim1vplrfWE5o/4mJhO/jqiK8aMO9LsEJzq+w372J3b2GB/dFoCUS6YRnNwQi+SIkMAsGqj6PXUXvGSY+aKjPK9MRJFo56UX2WVNVQ03P8Z6M+wlASTI4KAANdffkqMiSAk0Oi4Wl1bx5pd8oGyuzpSgO60fU1rtjwN2K5bvuOegFGsfgTU2h7zbD+rBe7vWqjCVeY9/ojZIThNem4Jz3631v58YHwEgxLaH2qpK5RSTBqQgMV2LSCntIr0vOZ3qggBO7ZtNTsEIZwiu7Dxf9yofgmm3/8JUFZa0v5K3WRRigGJjR8kf9ywx+XH9HUdLUAzgbMaFiilAoDTgB9aWf8L4Khmj6dsPzsKeLnr4QrRtnqr5t7PVlNt6xAUGRLA0Wmu/XQeFRrIyD6N/5TWZuRTb/XMVlAhhOiuLIfL72PTEk2MxP3SkhpbtNelZ9lbgkXXtHsPqNZaK6UeB55VShUCvwK3AHHA0wBKqYFAvNZ6hdY6H2gyXYBSapptX6udHL/oAaoqK6iormmyrKxW4x/gT11tnX3Z56t3sS7TuDncomDGkKQmwy211pO26V67ZnRyNHvyyymrqqGsuo4d2cUMT4pywp6FEMJzaK2bdEAaO6BnFaDhIUEkRIWTU1SGVcOerAJG9uttdlheq0NzwWutn1dKhWDc13k7sA44WWudblvlPuBKGjsmCS/iylkknKGiuqZFD8fgmgpmTxtj7xlbVlXLZ+v22X8+OjmmRYeg1nrSzp42ptvxBfr7ce7RQ3njp40ArN9f4LLL/sI7hYeF02ywECG8TmFZJTW2oY2CA/zpF9fLZUOXdYa/f4dKGadIS4qxF+HphwoYkZrgk9OQukOH/2pa63k03svZ/GdXAVcdZtt/Af/qVGTCbc4+/0KzQ+gWrTXL03PsQy6lxEY0mbvdHWYdMZAPl22loqaOqtp6mSFJNNEnJQX2tBwmpieyai0tFV7Ksfd77+hwaurq+WhFy/ub3T2UmDuP1TcukgA/C7X1VsqrasgpKqd3dLjbju9L3PexQXisTz9836uL0PS8Ug4UVdif3zRrPNv3HHBrDEEBfoxLibEPy7TpQCElFdUQ4vrhQVzNWwcBb+3WDXDe4NSdcTCz4zO1tDUgvhlxd5fWml2HClm9N4/s0kqKK2uoqbNiUYpvtx4iLMBCv9hw+kaHYZFWJI/n2AHJk4quqsqK9ldyEn8/C/16R7ProHGnYXpWgUe9Ft5EClBB+u6d7a/koSpr6li5J9f+fHhSFEOTY91egAIMSujFpoOFlFTWUltv5bPVuxh18ii3x+Fs3joIeGu3boDzBqfujLLysvZXsmlrQHwz4u4qrTX78svYcKCQN5Y3n0TPaAXNKjJG8NudW0p4kD8T+sfRL0beyD1VXb2VXIf533tHNx+Z0Tx1dXXtr+REaUkx9gJ0f24RNYOSCXTCJCc9jfnjJwjRDSv35lJdZ7TChQX5c2RqbDtbuI5FKY5IaTz+gjW7Ka927z9GIcxWWF7N15v289OOLArKqzu0TVl1HT9tz2Lx9kOUS89ij5RTVGYf4SMiJIiw4J479Vt0eAj9EyIBY/SVfTmFJkfknaQFVHit1bsOsSevsWVpSlpCk17vZugXG05EcAClVbWUVdXy0Zp0rpziPS1XQnRVTV09q/fmsflQIY6jQwf6+5ESHUq/2HDiwoMJDvCj3qoZNTiFdxavZ2d2sf1DZEZBOXe9+SOTB8TSqwcXOJ7oUEFjJ7o+sT27k6VSiuPHDGD+onWA0RlpcHKcuUF5IWkBFcy9+16zQ+i0mrp6Xvq+6YDzydFhJkZksCjFKIdxQV9ftpOaOs+9T1K4x5Bhw80OwaUyCsq4/c3FbDrYWHwqBaP6RPPSjacwfXAiqTHhhAb6Y1GKAD8LQ/vEMqFfHOcc2Z+hvSPt+zpUWMY3mw5QXOmMQdKEM2itOZjfeP9nUqznXH4HCI9wf0E8fWQ/+33LhWWVFJZVuj0GbycFqGDDut/NDqHT1uzLJ7/USPhgfz+O6h9vckSNBiZEEGy7Hyi7pJKFGzNMjkiYrbjIdy/Rfb/lAOe9sIjd2cX2Zb17hXDm2FQm9I8jsp1pcIP8/Zg8MIEZQxLxs00rVlFTx7ebD1BRI7eweIID+SWUVxkfCPz9LMRHmv9h31Ftrfs/rESEBNI3vvGD055DBW6PwdtJASr4/puvzA6hU3KKytju8GY3KS3eXvB5An+LhREOA9G/+st2rDI7Uo+WnZVldghOV1tv5R9fr+e295ZTZrvX2aIUkwbEM2tkMlHtFJ7NDYiLYObwPgTZcrmipo5FWw9SWy9XEMy2amfjvOeJ0RH4WTyrdKiuqjLluGmJjTMj7c0upN4q52pneNZZJEQ76uutrNq+3/48JTqM/rGe13N2aGIkoYHGLdbpeaUs3n7I5IiEcJ5DxRVcOf8n3ljeOIJG78hQTh3dl+FJUV0emDspMpQ/nzXZPk5oQXk1y3fnoLV8gDPTGocCNCnGsy6/m6l3dDihQQGAcVvYgbzidrYQjny6AK23WskprWR7VjHrMvNZn1nA12t3k11SSZ18UvFKm/ZlU1pp9KwNCfTn6LR4j5yFIsjfj1OPSLM/f3uF9w51JYSjpTuzOPeFRfZpbwGOG5bEs1efQFx498e9HTegN5MHJtifp+eVsmjD3m7vV3RNeVUNmzNy7M97egckR0qpJvPD75bL8J3ik73giytr2HqoiA/WfElFsyE91mYaY3f5WxSpMeEMS4okISLEjDA9xpnnXmB2CB1SWFrJNod/hFfMGE11J8ZXdLczJwzk01U7qbdqfnMYq1T0PH2SUyBnX/srerDaeivPLd7Cf3/ZZu9o5GdR/HHmKOZMHUKmE6+CDukdSW5pFTtzjJ7X839cz+mjU4gIDnDeQUSHrNqeaR9+KSo8mJAgz/sbuHMmpOYGJMawaW82YMwUlVNUTkxUlGnxeBOfKkCr6+pZl5HPtqxi2rtgU2fVpOeVkp5XSlJkCEeP6EdqapQ7wvQ4iYlJZofQLqtVs3J7pv3v2rtXCCeOG8CXv240Na7Die8VygnD+/DdZvcPiu8OnZ0hqaDeQllN09l8/AP8qatt2dHE02dZ6qzgYPfPiNXaTFBdnaM7s6CMOz9ayYb9jS088RHBzDt/EhPsHQCde5l80oB4ckurKKqsobq2nqW7spk1Mtkjr3j4sp83ptu/7xPjma2fZt6TGhYcSGJ0OFm2aUoXrdvNsP7JpsXjTXymAM0treKnHYdaDPwdGuhPYq8QwoOMXzU6MoLf0w9RWtXYMnqouJKbXv2Ba6cP5ebjRpg+lqS7vfTcMx4/FNP2/bn2YS4sSjFlYIJXTN132aRBPluAdnaGpLKaet77uen6s6eN8cpZljrLjNnGWpsJqrNzdFutmg/X7OHJ7zY0+d86OS2BJ86bSKwTLrm3xd/PwrTBvVm4wfjgmV1Syc6cEoY4DNkkXKumrp5lWxpb7lPiPfO1Lzf5SlhaUqy9AP1u7W5uPGMq/n6e0zHWU/lEAfrz5gy+3pSJY0fjpMgQbpw1gfSMQ00+Mc+eNobPf1lPQXk1Ww4VkZ5bisaYGu7lJdtYvjubf188hd69evZleU9ysKCUTXsbexGP6t+byBDvGKR6fL84hiVGsS2ryOxQhBfQWlNQWklecTmVNbXU12usWqPRLN+2n+zCMoID/QkLDsTfxR+Utxws5MEvfmfjgcYhpPwtittOGMXVU4dgsbj+A0JceDCjk6PZYIthzb58+sWGE+Qvb+7usHpHpn34pbDgQKLCfft9UWlNQVFRk2X1dfWtr+wgOa4XwQH+VNXWUVBaydJNezl27EAXRek7vL4A3XKoiJXLGlsXAvwsTBmYQP/YcEb3S2BPZsvhT5RSxIYHM31wIqP6RLNiTy7ZJUbr2sYDhVz08o88d+kURiRFt9hWuJfWmue++K3JPUjDUhKgzjsG/VVKcdnRg7j3s9VmhyI8UFF5FXuzC8kvqaCgtIKiskre/7llizDQZPQHMAqCyNAgrHV1jBucwtC+8SREhXf7EvW2zBzeXbyWn9bvbvKhPjUmnH+eP5HRyTFtb+wCY/rGcLC0mrySCuM2q8wCJg3wnHF/fdlPGxovv6fER/r87Q/VtXV8tGJHk2XnHd3+THZ+FgtpfWLYss/oo/DJr5ukAO0Ary5At2cVs9Khc0dUSCAnDO/TqRvVo8OCmDUyGb+QUF5fvBGr1mSXVHLpKz9x95kTOXFkik/di9aa0WOPMDuENn2xYisbbUmtgIlDU9zS8uJMp41OYd53G9hjdiDCNJFRUZBTgNaawrJK9ucWc6iglPd+Wt/lfZZX1VBeVcMHSzfzwVLjUntUeAjD+sYzNCWe5JgwKqpqCAkKaLdwKC6vYtmWvXyzejtrdja9ZSTAz8K104dy3fRhpoy36+9n4apjR/Pkgt8A2HaoiCG9exHdyXFGRefU1dezdFPjf62+cZ55+R0gIMD8K2KDkmLZui8HDfy+6wB7sgoYkOjeD2vexmsL0L35pSxPb+wRHR8RbAxi3IVLM0opTp8wmP1ZeSzenkVtvZXq2noe/Gg55VW1zJnY34mRe56TTjnN7BBalVVQyrMLfrU/H5oST0yEeb0duyoowI/zJ6ThffNNCWfQWlOmg1i/+yAZucX2S5ptiQgJxN/Pj5CgAAL8LCilUAriI0JIzymmqrqW8qqaVrv8FJVVsmJbBiu2Nc6+FRzgT3RECGHBgdTX1REX1YuqmlrKKmsoKKsg/VABe7JaHz7mmMGJ3H3qWPqbPPXi5KHJJPYKIaukEg2s3JPLSSOko4crrd19kJIKY8i72IgQYnt57v/eIBM6+TUXGhxIclwk+21jgX766ybuOPcYk6PybF5ZgOaVVfHLzmz784GJ0UxKjSawm/cF9YkK47TRKSzaesA+s8e8L1cTF6iYPa5ft/btyd567RUuv/pas8NoQmvNY+8vtg+jFRESxKj+iSZH1XUXHpXGPWYHIdxGa82ug/n8uG4XP67bxYH8klbX87MoYiJCSYgKJ7ZXKDERIVx2zMgWlwHBuBTYsLy+3kppZTXF5VX0jgxlb24JO/bnthh2DqCqto5DBcY83rsO5rcbu0Upjhs3kD8c3c9jbkNSthmWFqzPQGN0HN1fWGF2WD7tZ4fL75OHp3j05feK8nKzQwBgUHKsvQD9ZvV2bjjtaMKCzW+d9VReV4BW1dbz47ZD9nsCewUHcP/5U1m8ZrtT9h8VGsjpY1L4dvMBCiuMVoZ7Pl0F4LNFaE62500T+NmyzazZadzzZlGKScNSXN7pwpWSIj239UA4R0PRuXj9bn5cv4v9ua3PihLgZ6FPXCQpcZFcf9I4vlrb+Zsz/PwsRIWHEBUewiXHjCQmKgqrVbM/r5htmTls35/L5r2H2LY/j7oOTGXpZ1EMT01g2sgBzDxiMIkxEcSU5bS7nTtFhwUxNDGSbVnG67ouM982Q5LnFkbeqq6+vsn9n1OGp7Jhn+eOZWy1tt9RyB16R4XTN64X+/NKqKiu5cvftnLhjLFmh+WxvKoAtWrNkp1ZVNQYrZOBfhbjns8Q594LFBzgz8kj+/Lt5v1GEaqNIjQowI+TR/Z16rFESwfyi3n+i2X252dNHkZgoHyKFJ6ntLKaNTv289u2DH7bnkFOUestMcEBFnrHRJIaH0liTONc2qFOHNTbYlGkJkSRmhDFSeOHUFBUxDs/b6K0spqisiqqa+sYlhxLPYqQwADCQ4KICAkiJT6SgUmxhHpBS82YvjHszCmh3qrJL69m2Y6D9Bsn/5Od7bdtmRTZhr2LjwxjeEqcRxegnkIpxexJQ3l+odFo9cGS9Zw7bbRXN564klcVoB//toODRY2XXaYPTnTZcDzBtmLzt3357MkpRmu466OVRIcGMnFAQvs78CJh4Z4zt6/VqnnsvcVU2j5k9O8dzSXHjuGjZdtMjkz0dFarJiO3kK0ZOWzJyGHLvmx2HcyzX41pLiQogOkj+3P8uEH4leeyJrvl5XFXU0rRKzSYXqHGPXIX2lpLvVVooD/DEiPZfLAIgLd+2cKFY5K9rmOip/vW4YriiUcO6fRA710dzqirlPKcAu+4MQN496eNFJVXkV1Yxk8bdjPziMFmh+WRvKYA3ZVTwptLttifj06OJiUmzKXHDA7w4/GLp3P3Oz+zJ6+U2nort7y7jLeuOZahiVEuPbY73XjLbWaHYPfx0o2s230QMC4L3nvJCd2+t1eIjtBaU11TR3lVDWW2HublVTXs3J9DfmkVhwpKqW5l1iZH4cGBTBqeygnjBjFpaCpBgca/2IKiaNZkbz7stqJjRidHsz2rmDqrZm9uCd9s3s+po1PMDstnlFZWs3TTXvvzWRPaH4aoua4OZ9RVYeHhLtt3ZwUF+HP21FG89p0x9N7/Fq/jhHGDPPoeWrN4RQFaV2/lr5+uotZ2L1NsWBBHpMS65diRoUG8fPk0LnllMbmlVZRV13HDW0t557rjSI5ybQHsLsuWLmHKNPN76+3LLuTFhSvszy874UiGpSS0+CQtRFfV1NaRnlVIZm4R5VW19iKzvKqGz5Ztpqrm8AVma4alxDNpWCqThqUyIrV3q5fb8vPk8qWzBAf4MyIpyj44/XOLt3DSiGS5zOkkP63fTY2ttXJIchxpSbEe/z+4prra7BCaOGfqKN75cS01dfVs35/Lut0HOWKQjNrQnFcUoPN/3WGfjcOiYNqg3m695JIcHcZLl0/jild/oqy6jpzSKm54cylvX3ssUT4wFt1yDyhAq2vruP/N7+wtTAOTYrnqxAmmxiScp6aunqoaKzV19Sil8Lcop7cIWK2agrIKcovKyCkqY8+hXH7fdYDyysYis7beyie/dr0lMjo8hOGpCYzo15sRqQkMS0mwX94+nPy8vC4fU7Q0sk80W7OKqa23sievlO+3HOAUaQV1im9XN7ZcnjxhqImRdFxNjWcVoNERoZw8YShfrDCu2r71w+9SgLbC4wvQndnFPLu48Q1jXEos0WHuK/qsFgsZlVZCIyO579zJ3Pv+r9TWW0nPK+Wmd35l/lUz2h2cuaqygorqlmP/hQYFdmpeZl/23IJl7D5kDBET6O/H/ZfOJEAuvXud8upaft6Syf7sArZnFZNTWklOaZW946Ajf4vis/WZWOvrCfCzNHnkVdbTO8yfsKAAAiwWNBqrNjoilvrtpbSimtKqauNrRTW5xeXkFpdTb+3epBH+fhbCggMJCw4k3PZ11viBDElJIjEmwukdHkXXBAX4MSIpivX7jfFLX1m6nVmjpDNSdx0oLGddunELlEUpuXexGy4+dhwLf9uKVWtWbs9k875sRvbrbXZYHsWjC1ANPPjF79TVGzf5D0mKZlSye8elq6q1smBp49R4UwYmsGRHFhpYl1nAnz9eyVMXHI3fYVpkK6preHdJy1aXS44ZKQUosGRjOp/8usn+/NYzpzKwj3tusRDdV1VbT3peKbtySigor4Y1ezu0XZ1VU1Re1erP0vNKnRhhIz+LIik6nHpNkyIzLDiQi6cP55t1e1u0zE4eluLVHXd81fAkY0im6rp6th4qYtnu7PY3Eof14ZrGIcGOGpri0YPPe7rUhChmHjGI7343pgp//btV/PO6002OyrN4dAFaE5PG7xlGq5i/RXHHaRNYsWWvqTENiItgeGoCLy4yitLvtxzgn9+u5+5TxpkaV3dceuUc0459qKCEx95bbH9+zOgBnDVlpGnxiI7LKipj6a5s0nNLaKMjuJ2fReFvsdhbKOva26CLeoUGkRAVTkJUOL1CAjhYUEZ4SGORGRTgz/mTh7Y60HtESJDLOgqk9hvA6pzdLtl3TxUc4M/JY/uzYI3xur7yi3PGgu6pauqsfOxQgJ45eYSJ0XROaKhn9se44sTxfL92J1rD8q0ZbMvMYViKb42i0x0eXYBWpjTeA3jF5MH0j+/FisOs7y5nHTWYivIK3ly+C4A3l++iT2QYV0yRyxWdUVVTy19e+4bSSuP+nd7R4dx94XHSW9DD1dZbWb+/gLdX7GpRSFqUYkxqHJP6xTKmbwx9Y8JIiAihV3AAmVXw3s/GBzer1tTVWzlu/DAWLt9MXb2VWquVmjordfVWRqclEWito6y6FqtVY7FNR2lRiprgcCJCgogIDbKPZRnbK5T4yDCCAxvH1SwoKmr1yoPwHedMHMzCtenUWzW/7ZGOXt2xaOsB8suN/8UJUWFMGdHf3IB8QP/eMRw/dhA/rDNqhde+W80/rjnV5Kg8h0cXoDogBICEiGD+cOxw8rt3e5dT3XXyWLKKK/luywEA/vHtehIjQ5gwyfs+3bzzxnzm3n2vW4+ptebx939i5wGjc4a/n4UHLz+pQx06hHmySypZujOb0mZTPsaHBzMooRf948K58oQjSA1p7UNEY7FqUYpAfz9iI0KICm05lu/s8QPb2AcUhHtfjmXs6/xsR6J9iVFhnDIqhS83ZJgditf738rGFvozJo3wqlEFKio8YyrO1lx50gR7Afrr5r1s3HOI0QOSTI7KM3h0Adrgz7PGEhYUQH6lay7bdYXFonj83Inkli5hbWY+WsOfP17JvxKS5OTqgPd/Xs+itTvtz28/Z7pXz/Xu67TWbMsqZuWeXByzMCEimPH94ujdK8Spx2vo/NdceKBndEzzpI6FrQ36DZ0b+Lu136esVhMe6EeMX8u/Q0G9hbKapvtv+Wq4x7XTh0oB2k07s4tZs89oDPCzWDjjaO+5/O7p0hJjOOGIQfyw1ihCn/9iOc/ferZc6cMLCtBJA+I9tndjcIAfz14yhUteWcy+/DKq66zcPf9rXrj1HFIToswOz2Mt37KP579Ybn8+++gRnDlZ7vv0VHX1Vlak57I9u3Fu8wA/C9efdARlRcUu+UfavPNfg4tmjMEThpz2pI6FrQ36DZ0b+Lu13ye4poKLZowhppWW6LKaevvtFA1mTxvT4eM505DekcwYksjPO7JMOb4vcGz9nD56AHGRnnlPpbe6/tRJ/Lwhnbp6Kxv3ZvHLpj0cMzrN7LBM59lt7NrKX04b59GfFKLDgnjp8mnE2IaGKi6v4k///ZLC0op2tvQck904BujWjGzue/NbrNpoRxvVvzd/PGe6244vOqe23sojHy1tUnzGhQdx5rhUThjd36Nz05PExsWZHYJPu2aad4xX6Ynyy6r4dO1e+/OzvbATaGCgZw+PlhwbyTlTR9mfv/DlCurqXTc1qbfw6AI0MHcngxMizQ6jXakx4Tx3yRT7eKAH80u469WvqKpx/9zPXXHk+AkUFBW1eFRVdryIrqqsaHcf+3OLufOVr+yzzSTFRPDoVafIVJseqrbeyqKtB9mwr7Fzx4C4cGaN7Et4UMBhthTNxcbFmx2CTxvfL85ts+P5mrdX7KK6zrjNYkRSFEd64YDpgUGeXYACXHniBMKDjfvdM3OL+GyZdJDscAGqlLpOKbVTKVWplFqulJrczvpTlFKLlVJFSqmDSqk3lVKdGoU15OC6zqxuqrEpsfzzvElYbC1CWzNy+Ovr31LbifuwzLJp8xbeXbK5xaO1e9za0nAJr6195JdU8Kf/fklRWSVgDJcz7/rTZZw5D1VXb+XHbQfJLqm0LxuXEsMxgxO9qnOCp0jftcvsEHyaUoprpksraGeVV9c2ufx+7fRhXnlVo7yszOwQ2hUZFsxlJxxpf/7K1yvJK2t9HOSeokPvJEqpK4AXgbeBc4Ei4Ful1IA21h8O/ACUAhcDfwKm2rbpcNOJpc67/jgnDO/DbWdPsz//bVsGD779fbdnZ3G1ujrXttTml1Rw2wufsz/PuIwb6O/HP645ldQE904qIDrGatX8uP0Qh4obi8/x/WIZlxLrlW9OnsDVOSbg2CHS+bOzPli9h5Iq49xMjQnnxBHe1/oJoLVnv8c2uGDGWPrGG1d1y6pqeOq7jSZHZK52OyEp4x3nIeBlrfWDtmXfA9uB24H/a2WzW4BDwLla61rbNjuBlcCJwFdOid4DnTttNEVllbz23WoAftqQjp+CuKjwHvnmXVhWyd/e+Zq92YWAMSD5A5efKCMFeCitNSv25HCwqPHWiUuPGUlAnVl9nIXoGMthZqMTLVXX1vPGssbOa9dMG3LYGf1E9wX6+3H72dOZ+/KXAHy2bh/TRvVnVErTe8TbGn3C13SkF/wgoB+woGGB1rpWKbUQmNXGNpuBLQ3Fp03DNBWttpr6kjknH0V5VQ0fLDF6if6wPp20pBiOGtLXI4vQoKBgwPmdpiqqa7n3zR/IzCsBjOLz/stOlN5/HuyTlTvZkV1ifz62bwznHj2s1R7pouNclWNCdNX/Vu4mp9S4yhgfEcyZ4/qZHFHXWSze049g0rBUjh2Txk8b0gF45JPlnDEmtUnx39boE76mI5fgG8byaH4TUzowUCnV4i+vtX5ea/1cs8Vn2L5u61yI3kcpxa1nTuX0ScPty9IPFfDbtkx7729P0m+A8z8TFJVV8v3vO1sUnyeMG+T0Ywnn+GHrQV75sfGSUFp8BONSYkyMyHe4IseE6KqSyhpeWtL4Vnzd9GFe3Rk0NMy7ho265cypBAca7X9FFTVs2F9gckTm6EgLaC/b19Jmy0sxCtgwoITDUEqlAE8Cq4EfOxmjV1JKcef5M6itr+fb1cZljr3ZhVi15uhhqR51uSg765Bz91dYxtJNe6itNy4h+Fks3H/ZTCk+PdievFL+/PFK+yDzCRHBTB2Y0KUW+7YGkXflRfzWBlLvzEDsrubsHBOiO15dup3iSiNfUqLDuGCCd1+Vqq7yrv4iidERXH/qJP792a8AbDhQQL/YcPtwjj1FRwrQhneg5k13DcsPe6OCrfj8AaNYvUjrlk2ASqnrgeubL09J9s4bohv4WSz85aLj0fX1fLfW6GmYkVNEbV09U0Z4zuWO4lZmUekKrTW7DuazdtdBe0tvSKA/j149i4lDU51yDNE1h8uxypo6bn9/BRW24bHCg/w5flgSfpau9XZvaxB5Vw5U3tpA6p0ZiN3VnJVjwjN503tYTkklb61ovKD5fyeMJNDfu0e2qK31vnvUz502mkWrtrHlgDGT4q+7sjltTIp9JJ2eoCNnXcMI1BHNlodjFJ9tTsKqlBoFLMNoRT1Ra727tfW01i9rrSc0f8THeP/lPz+LhZtOn8igPo1j1B0qKOXHdbvJ96LB6ttTWVPL8q0ZrNl5oEnx+fhVJ0rx6QEOl2OPfrWOHbaB5gP8LBw3rA/BAR4/SZoQHsOb3sOe+n4jVbXG1YHhSVGcMirF5Ih6Jj+LhdtPO9JecOaXV7M+s2ddiu/Iu0zDhN1pNL0PNA3Y3lqLJoBSahLwNcbl+eO11jtbW88bdXaeaotSjB+cTKC/H1sycgCjd/idr37HP687ncHJnjlLSltzTDef73r7/lweePNb+/2eANHhIUwb1Z8Bia0PteRJc2n3ZNWxA/nk97325zedNI7i4uZ32wghfMHy3dksWJ9hfz73xNEedTtYT5MS24sjUmNYsy8fgA37C+gT1XPe/zpagGYCZwHfAdjG8jwNWNjaBkqp/hjFZzZwgtb6oBNi9RhdmadaKcWYtCTCggNZvWM/GsgrqeCGZz7m9nOmc/qk4ab1kE8bOJjVOdtbLG9rjumG+a5r6up57btVvPvjWuqtjZ9D0pJiGD8oGb/DDFjuSXNp92QVqZPs388em8qssf15f0nPHpvOFdrKMSFcraDeQllNPTV19dy34Hf78hNGpjBlUMu5YTz9furWhIW19c7rfh1tuGkwsk80BworyCqpRANLdmRxVXUdGa2UZ742PFO7BajWWiulHgeeVUoVAr9ijPMZBzwNoJQaCMRrrVfYNnsG47L7zUCqUsrxGuw+rXWPvSN/YJ9YQoMD+XXzXurqrdTU1fOPD35iXfpB5p47g1ATpjis6sIN3L/vPMBTnyyxj+8JRk/3CUP6MiDR8y47iTb4GedbWlwE959xJHn10hriCl3JMSGcoaymnvd+3sDajHwOFhp3zAX4WZhz/OhW1/f0+6lb40mTvbTXcNOcRSmmD+7N5+szqKmzUl5TxzML1zA4JrhFo5SvDc/UoRu9tNbPK6VCgNswBp9fB5ystU63rXIfcCXGuPUBwKmAH/BuK7u7E6NHfI+VFBPBSeMHs3lvFvtyjHvvvl29g/W7D3HXBTPcfs/kwQOZHV63tKKaR99fwm/b9zdZPiI1ngGJsUSE9qxefL7A36L4x3kTCQ30h0rPGybMF3Qmx4RwttzSKjYcaLy/cHy/WGLDQ0yMyLmqKr27P0VYUABTB/Zm8XajbW7Z9v2ogQkM7h1pcmSu1eGub1rreVrrVK11qNZ6itZ6ucPPrtJaK9v3tVrrAK21auPRo4vPBr1Cg3nympM59ahh9mVZhaXc8dKXPPq/HygsqzzM1u5XWlHNym2ZfLVqW5PiMzjQn9vPns7fr5wpxaeXumHGcEb2kWlRhfBFZVU1/LzjEA29NRIighnq44WNN+oXG86Q3r3sz1fsyaWgvNrEiFxPurqaKCjAn79cfDzjByfzzGdLKakwTravV23npw3pXHzsOC6aMZbQ4EBT4tNak1dczo4DeWTmFLUYh+vkCUO44dSjSYgKb/WeF+H5/MrzuP6YYe2vKITwOlpr/v3NWsqqjSHWAvwsTB+c6JEz8gk4qn882SWVFFfWUm/VLN5+iNPHpBDkxZMEHI4UoB7g5AlDOWpoCs98upQf1hkDDVRW1zL/21V8snQj50wbzezJI4jr5ZrZHnonJkJOYz+xqppaMnKKmbt5Dxm5LecYGJmawB/Pmc7w1JY3sAvvErZnKQF+N5odhs9rnmNmaauDhKd3MhFd8+bynSzZ2njFaurABCKC3d/PwNWCgoPNDsEpAvwsHDc0iW82H6Cqtp7Sqlp+2ZnNCcOSfPJDgxSgHiImIpQHrziJkycM4YUvV7Any7hfp6i8ivnfruKN79dw7Jg0Thw/hKMG9yUo0Hl/ul6RUZRW7CGrsJTM3GJyi8patHaCMXvDiH4J/HH2JGKiopx2fGEev6ri9lcS3RYZFQ2YX4C21UHC0zuZiM77fssBnvi2cbSWIb170T+u+XDeviEgwJyrhK4QFRrELadM4MkFvwGwv7Cc9fsLGJcS286W3kcKUA8zZUR/Jg1L5bs1O3j1m1VkFRpjMtZbrfywbhc/rNtFcKA/E4emMG5gH0ak9mZwchxBnRg4vLK6ll0H89hxII9tmTks35xOUUVtq+v6WSz06x3FoD6xxETI8EhCdMWObVvNDkF4gIYhkZrr7PA6re3HcR/rM/O566Pfmtz3ObF/fJP12xrPul55X2t4WelhZwP3CG1dfWjNlGF9Gbl6O5sPGuuvyywg1gen6ZQC1ImsFkunLm+1dkI2jBV2ylHDmHnEYJZsTOeTXzexPr1x5KqqmjqWbNzDko17AGP4o4SocHpHR9A7KpyQoAD8lPGfp6aunqqaOsqraigoqyS3qKJDMzDFR4Zx1qSh5FfUENiN+086e8mvtfXl8qAQwhc0DInUXGeH12ltPw37WJ+Zz/VvLaW6zigu+0SHccyg3vg3G5e5rfGsTz9mXIfjEB3X2tWHCyYOanP98f3iyC+rJqvE6JD8844szsoqJHWA7wxzKAWoE1XVWvng15aDq7d1eau1E9JxrLAAfz9OOGIwJxwxmF0H81j0+06WbNpDRk5Rk23qrZpDBaUcKuj6DDb+fhYSIsNIiu1F37hIQoICOPnIga1eruuMzl7ya219uTwohBDtW7Unlz+88ysVNUano6jQQB6+YCq/2BorhPewKMWxQxP5ckMmZdV11Fk1f/twGSNuOJ6kSN+4GikFqJcY1CeOQX3iuPH0yezLLmT1jv1sychmS0YOmblFndqXn0WREh/FkOQ4BveNJ7SunH1lWqZkE8JFwsPCMWYlFsI1Fm3cx3+++d3e8hkdGsirVx1DaKRv3vfpyN/fN0uZ4AB/Zg7vw8KN+6mtt1JQVsUf3v6Vt685lnAf6Ezmm381H9evdzT9ekdzLsZMFlU1tWQXlpFdVEpucTlVNXUUlZaxetch/C0W/P2MxxkThzAkJYm4yLAml2MKiorIbGVaTCGEc/RJSYE9kmPC+eqtmlV7c9m2bKd9WXxEMPOvOoaB8b3I6AGTS/jy9M1RoUEcPyyJ77YcQGvYkV3M7e+v4LlLpxLo3+Gh3D2SFKA+IDgwwF6UNigoKqKytukN5qP6JRAT1fLT8MFMmaVFCFeSHBOukFVcyfL0HIorG+duHxjfi+cunUJqjOfMj+5q3j4TUnuSIkOZOrA3S3dlA/Dr7mz+/PFKnjx/En5efOVSClBBWXmZ2SEI4dMkx4QzZeaXsGRnFum5Te/7P3lkXx45awJhQT3rrb2urs7sEFxuUEIv+iXG8M5SY0SNbzfvJyTAj0fOmuC1t8/1rLNUCCGE8EI1dfX8uiubT9fu5YetB5uM1exvUVw3cxyzjxxAvlWR73DZvablroSXuuSYkZRV1fD56t0AfLZuH5Va8ZfZE4j3wltCpQAVQniNzg515gw9eWiwtsaKlKLG9cqqatl8sJCNBwrYsL+AlXtyKalqOV5zv5gwJg6I58Rxaby/pOWwSrOnjXFHuMINqus0UQGKwQm92JljdGr8dv1eamrr+df5RxHg5133hEoB6mHaGjcz0KKosba8mbyt5Z15gxwybDirczreQaInvyELc3V2qDNncMbQYJ3NMU/R1liRUtQ4l1VrtmUVsy4jnw0HCti4v5D0vBL7QPKt6RsdxpjkaBJ6hbgvUA8VHtELyDI7DLdQSjF5YAJWrdltuwVj8ZZM/vh+HU+eN4kQJ86S6GreE2kPcbhxMzu7vKOKiwq7HaOM1SlE2zqbY8L3WbVmf2E5e/JK+XTtviYdidqSHBXKKaNTmDA0lRVb9ro+SC9RW9uz2uQtSjFtUG8C/CxsyzKmU1687RBXvfYzz106lbjwYJMj7BgpQAXZWT3jk6MQZpEcEw1q661sPljI9qxiKmvbvnJkUTC4dyRjkmMY0zeG0X1jGJzQC6UUGZWaFW6M2dNVV1WZHYLbKaWYNCAefz8Lmw4YH3A3Hijkopd/5NlLpjAsMcrcADtAClAhhBDCxeqtmoVr0/n4971UtVJ4RocGctSAeMb2jWV0cjQj+kQT6kWXU4X7KaWY0C+OGSP78cL367BqOFhUwcUv/8i9px/BuUcOMDvEw5KzWwghhHCh7JJK7vroN1btzWuyPCTAj0EJvbjx5COZkRrltcPpCHOdMX4goxLCmPvBb1TU1FFdZ+W+z9bw+7587j1tnMfeF+qZUQm36pOcAjn7zA5DCJ8lOdZzLdlxiHs+WUVhReN9iqGB/hyZGktafAQWpRiQFMP+aivQtNdReKAfMX4tRyEQTfnyTEgdYbVYGJCSyDNXHccjn/zGvjyjh/yna/ey6UABj587keFJUQAU1Fsoq2nZAm/GuSYFqCA42DtuWBbCW0mO9UzVcUO46Z1faRioxKJgVHIMY5Kjm0yH3NZoAxfNGENMiLSKtsfP4l3DDzmb4/kzbWA8/qqxh/zOnBIufPkHbjp2BNdOG0pZjea9nz3jXOvZfzUBQPrune2vJIToMsmxnqmi/2R78ZkQEczjlxzDkamxTYpP0X3lMtOYXYCfhWmDenPbKUcSHOAHQF295t8/bObSVxaTmV9icoSNJAuEEMJJGsbIbf6IiY01OzRholF9ovn4DzMZkxpvdiiiB1BKccq4AXzyh5mM7RtjX77xQCE3z/+BzQcLsR5ukFk3kUvwPUhbg9zHxMZCTq77AxLCx7Q1ju+skYkmRCM8wbRBvXn6wsmEBflTXmn+m77oOfrHRfDWNcfy2q87+M/izdTVa2rqrKzam8fevDKmDkogKjTItPikAO1B2npzPGNcqgnRCNFzBAQEmh2CMEFA4T6eu/Qcr5si0dtIfrXN38/CdccM45ghidz9ySq22wauzy2rYsH6TMalxDCqT7QpsUlWCIKkg4QQLiU51jOFpS+R4tMNJL/aNzQxivevP4HLp4+gYbQvq9b8npHPlxsz2Z1d5PaYpAVUUFFebnYIQvg0ybGeSenuDWtjtVjIqGy5j5418WT7JL9aav3cUZw/bTglJWX8ujubvLJqAArKq/m/13/kuulDuXHGcAL9/dwSoxSgAqu17enghBDdJzkmuqKt4ZlmTxtjQjSeS/KrpcOdO9FhQZw6OoUtB4tYm5lPvVVTb9W8+PM2vt9ykEfOGs/YFNd3nJRrA0IIIYQQPYhFKUYlRzN7bCoJEY23MOzOLeHSVxbzxDfrqaypc20MLt278ApKyWkghCtJjgnhOpJfXRcZEsgpo/py04ljCQk0Lr1bNby+bCdnP7+IVXtcN0KO/NUEYeHhZocghE+THBPCdSS/ukcpxewJg1hw80lMTkuwL88oKOPK137m4S/XUl5d6/TjSgEqqKmuNjsEIXya5JgQriP55RzJ0WG8cuV0Hj5zPBHBAfbl/1u5m9nPfs/SnVlOPZ4UoIKaGkleIVxJckwI15H8ch6lFOeOH8CCW07iuKFJ9uWHiiu4/q2lTj2WFKBCCCGEED2YMWyTtj+qA4K56+zJ/O2co4kKdc1A/x0ehkkpdR1wF9AXWAfcobVefpj1RwHPAJOAAuA54AmtPWACUiGEEEIIAbQ9bNNFM8Ywc1Acj361jm827XfqMTvUAqqUugJ4EXgbOBcoAr5VSg1oY/0EYBGggQuAl4FHgbndD1k4W2homNkhCOHTJMeEcB3JL9eKDQ/mqQuO5t8XT3bqftttAVVKKeAh4GWt9YO2Zd8D24Hbgf9rZbObbfuerbWuAL5SSgUB9yilntFaO787lRBCCCGEcImZw5Odur+OtIAOAvoBCxoW2ArIhcCsNraZCfxgKz4bfAbEAEd1KVLhMhUVMo2ZEK4kOSaE60h+eaeOFKBDbF93NVueDgxUSrU2aeiQNtZ33J8QQgghhOiBOlKA9rJ9LW22vNS2fWs3X/RqY33H/QkhhBBCiB5ItdcpXSl1CfAOkKi1znZYfh1G56IIrXVZs21qgPu01v9wWOYP1AK3aq2fbbb+9cD1rRx+FLCpU7+R+eKAPLOD6CSJ2T2Gaq0jzDiw5JjpJGb3MCXHJL9M540xg3fG7bQc68gwTMW2rxFAtsPycMAKtHbzRbFtfUcRDj9rQmv9MkYx24RSarXWekIHYvQYErN7eGvMZh1bcsxcErN7mJVjkl/m8saYwTvjdmaOdeQS/E7b17Rmy9OA7W2M67mzjfXB6D0vhBBCCCF6qI4WoJnAWQ0LlFIBwGnAD21s8wMwUynleH/oWUA+xiD2QgghhBCih2r3ErzWWiulHgeeVUoVAr8Ct2Dcu/A0gFJqIBCvtV5h2+x54FaM8T//CYwF7gHu1lrXOP/XEEIIIYQQ3qJDMyFprZ8H7gQuBz4CooCTtdYNQyvdByx3WP8Qxlig/rb1rwf+qrV+0mmRCyGEEEIIr9ThueC11vOAeW387CrgqmbLVgNTuxEbtHJTtxeQmN1DYnYOT4ypPRKze0jM3edp8XSExOw+3hi302JudxgmIYQQQgghnKlDl+CFEEIIIYRwFilAhRBCCCGEW0kB2oMppa5TSu1USlUqpZYrpSabHZMQvkRyTAjXkhzzXl5TgCql1imlzrJ976eUukMptVUpVa6U2qKUukUppUwOswnHmJstD7LF/rr7o7LHcAXwIvA2cC5QBHyrlBpgVkyd5cmvb1s8+Tz25Nja4snngOSYOTz1PPbUuNrjyeeA5Jg5nHYua609/gEkAjVAL9vzB4Aq4K/ACbbndcBdZsfaVszNfvZ3QAOvmxSbAvYCLzgsCwDSgX+b/dp5++vb0Zg96Tz25Ni88RyQHPOMmD3lPPbUuLz5HJAc84yYu3Mum/7LdPAXvgJYYvveApQADzdb5zkgx+xYW4u52fIjgDIg18TEHWw7sU9ptvw/wA6zXztvf307ErOnnceeHJs3ngOSY+bH7EnnsafG5c3ngOSY+TF391w29RK8UupdpZQ+zONb26onA9/Zvo8E3gQ+aba77UB8s+k/PSXmhm39gfnAP4EDroyzHUNsX3c1W54ODFRK+bk5HjtvfH09+Tz25NicHHPDtpJj7fDG19dTz2NPjctFcTdsKznWDm98fU05l02upEcBG4ElwNG2x8kYn2oeA/pjNLPnAEe1s6/vgUxPjhmjiXozEAisw7xPjhfb4k1stvxa2/IWlwO84Zww6/X15PPYk2PzpXOgld9BcsxDYnbleeypcfniOdDK7yA55iExd/VcNuWP0yzQQuAxh+czbL/wVNvz8UAeYDnMPhpOuFs9NWZgGFABTLY9NzNxL7HF27vZ8utsy8O97Zww+/X15PPYk2PzpXOgWfySYx4Qcyv7cPp57Klx+eI50Cx+yTEPiLmVfXT4XO7wVJzdoZSy0LTHvdZa1yul+mHMK7/e4WfjbMFvtD0/GViktba2se9LMXrBfQQ864kx2/b1KvCq1nq5s2LshmLb1wgg22F5OGAFyt0ekY03vr6ecB5LjkmOdZQ3vr5mn8femF/Ojtvsc6AVkmNOZMa57K57QO8Hah0eu23Lx9q+bnBY9whgr9a6xPb8ZOBbWqGUuh14C/gSuFTbym8PjPlWoB9wv1LK33aPh+1XUG75ENDMTtvXtGbL04DtTn4dO8sbX19POI8lxyTHOsobX1+zz2NvzC9nx232OdCc5Jhzuf9cdlOzbh9ggsNjtG35vUAl4Oew7nrgU9v3ERjd/ZNb2WfDEAVvAP6eHDPwky3Wth793fF3cIhHARnA8w7LGoav+I87Y2klNq97fT3hPJYckxxz5flq9utr9nnsjfnl7LjNPgda+d0kx0yOubvnsml/IFvQHwKrm508NcADtudnApta2e422y/7L0B5eszA0Gb/BCZg9BL7wvZ9oAmv/U0YlykeBU4FvsIYTiHN284Js19fTz6PPTk2XzoH2vg9JMd8/Dz21Lh88Rxo4/eQHDMxZtvyLp/LZjSbOxoL/OLwPBHjly6wPW/R3KuUSgL+gXFPwnvAJNV0wP3VWus6VwVMF2LWWm9vvhOlVCWQr7Ve7aI4D0tr/bxSKgTj5Lkd42bnk7XW6WbE48AbX19PPo89Oba2eOM50ILkmFN56nnsqXG1xxvPgRYkx5zK/eeyiZ8QQoF64DaHZcG2X6QSuBzjfpWTm213FYdvpo7ztJjb2Nc6PGyAWbMf3vj6evJ57Mmx+dI54E0Pb3x9PfU89tS4fPEc8KaHN76+Zp3LyrYTIYQQQggh3MLUmZCEEEIIIUTPIwWoEEIIIYRwKylAhRBCCCGEW0kBKoQQQggh3EoKUCGEEEII4VZSgAohhBBCCLeSAlQIIYQQQriVFKBCCCGEEMKtpAAVQgghhBBuJQWoEEIIIYRwKylAhRBCCCGEW0kBKoQQQggh3EoKUCGEEEII4VZSgAohhBBCCLeSAlQIIYQQQriVFKBCCCGEEMKtpAAVQgghhBBuJQWoEEIIIYRwKylAhRBCCCGEW0kBKoQQQggh3EoKUCGEEEII4VZSgAohhBBCCLeSAtTJlFKvK6V0O48HzI7TE9lemz+ZHYfwbJJjXSc5Jtoj+dV1kl+d4292AD7oYeBFh+dvAjttyxvsd2tEQvgWyTEhXEfyS7iFFKBOprXeDexueK6UqgBytdYrzItKCN8hOSaE60h+CXeRS/AmUUr9pJT6stmyPyqltMPzvUqpx5VSK5RSRbafP6CUWq2UulgptUMpVaWUWqWUmtJsX+OVUj8opSqUUrlKqf8opUKbHf9lpdS3SqkSpdSTbcQZoJR60HasaqVUoVLqE6VUSrM471JKvaCUKrDt7w2lVITDOsFKqX8rpXJsP39FKfV3pdTew7xGCUqpN237LFNKLVBKDejUCy16LMkxyTHhOpJfkl/dJQWo55sLfA1cZvsKMAR4CHgAOBcIAT5USvkDKKVGAEsADVwA/Bm4EPig2b6vBtKB84D32zj+08CtwOPAScBfgROAfzVb7y9ANHARcC9wse1rg/nAVcCDwCXAIOCOtn5ppVQIsBiYZjv+5UAisEQpFd3WdkJ0geSY5JhwHckvya/Waa3l4cIHsA54vZXlPwFfNlv2R+NPYn++F9jUbJ0HMJJyosOy2bZl423P/4eRlEEO60y3rXOMw/ELgIB24n8fmNNs2TNAXvM4AeWw7BNgo+37IbZjX+Xw81AgG9jrsEwDf7J9fwNQBwxz+HkvoBC43+y/qzw85yE5JjkmD9c9JL8kv1z1kBZQz7ellWV1wGqH5w03hIfZvh4HfAfUK6X8bZ8qlwMlGJ/8GuzSWtce7uBa6wu11vOVUn2UUscrpW7G+EQX1GzVldqWYQ4xNcQzw/b1M4f9VgALD3Po4zBufN/l8DtUAL80+x2E6C7JMckx4TqSX5JfrZJOSJ4vp5Vl1Vprq8Pzhu8bPlDEYnz6uqGVbZPa2XcTtvtyXgDGAMXAWqASUM1WrWj23OoQTxxQq7UuarZO9mEOHQsMA1r757KzvbiF6ATJsZYkx4SzSH61JPmFFKBm0rS8BzfcSfsuBj7HSLrm8jq6E6VUJPAlsBQ4V2u9y7b8CWBcJ+I5CAQopaKaJXD8YbYpBtYD17bys+pOHFv0XJJjkmPCdSS/JL+6RQpQ85QA/Zotm+6kfS/F+OS1puGSglKqN/AOxo3Xezu4n2EYN2X/yyFxLcCJtPz0eDi/YnyanI0xphxKqUBgFsalmLZ+hxMx7q/Js22jgLcx7tXZ2Inji55JckxyTLiO5JfkV7dIAWqer4EXlDGjxM8YvfjGO2nfDwPLgA+UUvOBYOA+IAXj8kNHbQNKgfuUUn4YPRVvBsYCWimlmt0z0yqt9S6l1DvAv5VSYcA+4P8wLqXsa2Oz+bZ1vldKPYZxs/n1GD0mz+jE7yB6LskxyTHhOpJfkl/dIp2QzPMKxie5W4AFGL3j/uiMHWut1wDHY1we+Bh4FTgAHKu1PtCJ/RRjJEu0LcbnMC5/nI9x7kzqRFh/AN4DHrV93Qd8CJS1cewS4BiMfyAvYlyO6QecqbX+qhPHFT2X5JjkmHAdyS/Jr25RHSj+hegWpVQcxvhrX2itSx2WLwOytNbnmBacED5AckwI15H8cg25BC/coRJ4HjhfKfUixj0z5wNHY9wjI4ToHskxIVxH8ssF5BK8cDmtdTnGp8dwjEsXCzDuwTlDa/2DmbEJ4Qskx4RwHckv15BL8EIIIYQQwq2kBVQIIYQQQriVFKBCCCGEEMKtpAAVQgghhBBuJQWoEEIIIYRwKylAhRBCCCGEW0kBKoQQQggh3Or/AY3xy6NbN5FbAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import seaborn as sns\n", "\n", "fig, ax = plt.subplots(1, 3, figsize=[11, 5], sharey=True)\n", "\n", "for i, dat in enumerate([dat_saz, dat_pfz, dat_miz]):\n", "\n", " window=40\n", " \n", " ml_t_roll = rolling_mean(dat.ml_t_smooth, window=window)[window:-window]\n", " ml_s_roll = rolling_mean(dat.ml_s_smooth, window=window)[window:-window]\n", "\n", " T_anom = rolling_mean(dat.ml_t_smooth, 4)-ml_t_roll\n", " S_anom = rolling_mean(dat.ml_s_smooth, 4)-ml_s_roll\n", " \n", " alpha = gsw.alpha(rolling_mean(dat.ml_s_smooth, 4), rolling_mean(dat.ml_t_smooth, 4), 0)\n", " beta = gsw.beta (rolling_mean(dat.ml_s_smooth, 4), rolling_mean(dat.ml_t_smooth, 4), 0)\n", " \n", " dT = T_anom.diff(dim='time')*alpha\n", " dS = S_anom.diff(dim='time')*beta\n", " \n", " R = dT/dS\n", " \n", " Tu = [math.atan(r) for r in R]\n", " \n", " Tu = np.array(Tu)\n", " \n", "# Tu, R, p = gsw.Turner_Rsubrho(dS, dT, p=0, axis=0)\n", " \n", "# # ax[i].hist(R, bins=np.arange(-7, 7.2, 0.5), facecolor=blue, zorder=10, alpha=0.6, edgecolor='0.95')\n", " ax[i].hist(Tu, bins=np.arange(-1.5, 1.6, 0.1), facecolor=blue, zorder=10, alpha=0.55, edgecolor='0.95', density=True)\n", " ax[i].set_xlim(-math.pi/2, math.pi/2)\n", " \n", " ax[i].xaxis.set_ticks([-math.pi/2, -math.pi/4, 0, math.pi/4, math.pi/2])\n", " ax[i].xaxis.set_ticklabels(['-$\\pi$/2', '-$\\pi$/4', 0, '$\\pi$/4', '$\\pi$/2'])\n", "# ax[i].grid(axis='x', zorder=0, c='0.25', ls='-')\n", " ax[i].plot([0,0],[0, 1.2], c='k', lw=1, alpha=0.5, ls='--')\n", " ax[i].set_ylim(0, 1.2)\n", " \n", " ax[i].fill_between(x=[-math.pi/2, -math.pi/4], y1=0, y2=1.2, facecolor=red, alpha=0.15)\n", " ax[i].fill_between(x=[ math.pi/2, math.pi/4], y1=0, y2=1.2, facecolor=red, alpha=0.15)\n", " ax[i].fill_between(x=[-math.pi/4, math.pi/4], y1=0, y2=1.2, facecolor=lightblue, alpha=0.15)\n", " \n", " sns.kdeplot(Tu, ax=ax[i], zorder=10, c=blue, lw=3)\n", " ax[i].set_ylabel('')\n", " ax[i].set_xlabel('Turner angle', labelpad=20)\n", " \n", "ax[1].text(-1.5, 1.0, 'Anti- \\n compensated', fontsize=13)\n", "ax[1].text( 0.1, 1.0, 'Compensated', fontsize=13)\n", "\n", "ax[1].text(-0.5, 1.25, 'Sal dom.', fontsize=16, color=blue)\n", "ax[1].text(-2, 1.25, 'Temp dom.', fontsize=16, color=red)\n", "ax[1].text( 0.7, 1.25, 'Temp dom.', fontsize=16, color=red)\n", "\n", "ax[0].text(-1.5, 1.13, 'a', fontweight='bold')\n", "ax[1].text(-1.5, 1.13, 'b', fontweight='bold')\n", "ax[2].text(-1.5, 1.13, 'c', fontweight='bold')\n", " \n", "# ax[0].set_ylabel('Probability')\n", "\n", "plt.savefig('../figs_submission2/fig9.png', dpi=300, bbox_inches='tight')" ] }, { "cell_type": "code", "execution_count": null, "id": "592146be-ce22-4924-9b45-b882dee7e1ee", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "id": "8807f23a-5328-4bab-b47f-1a521e615567", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "duplessis2021_JGR", "language": "python", "name": "duplessis2021_jgr" }, "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.7" } }, "nbformat": 4, "nbformat_minor": 5 }