{ "cells": [ { "cell_type": "markdown", "id": "45a684da-23b5-44e5-adb7-176048022777", "metadata": {}, "source": [ "## Plot the Yearly Flux Cross-sections for Final Figure" ] }, { "cell_type": "code", "execution_count": 1, "id": "919196bb-adf5-4288-bdc7-5f9e4112eb4b", "metadata": {}, "outputs": [], "source": [ "import cmocean.cm as cm\n", "import glob\n", "import matplotlib.pyplot as plt\n", "import matplotlib as mpl\n", "from matplotlib import gridspec\n", "from matplotlib.patches import Rectangle\n", "from matplotlib.collections import PatchCollection\n", "from matplotlib.colors import LinearSegmentedColormap\n", "import numpy as np\n", "import scipy.interpolate as interp\n", "import xarray as xr" ] }, { "cell_type": "code", "execution_count": 2, "id": "8188a53e-5d77-4f6d-a4ee-afd31157811f", "metadata": {}, "outputs": [], "source": [ "plt.rcParams['font.size'] = 36\n", "\n", "# get colormap\n", "ncolors = 256\n", "oldmaps = ['Purples', 'Blues', 'Greens', 'Oranges']\n", "newmaps = ['Purples_alpha', 'Blues_alpha', 'Greens_alpha', 'Oranges_alpha']\n", "\n", "for oldmap, newmap in zip(oldmaps, newmaps):\n", " color_array = plt.get_cmap(oldmap)(range(ncolors))\n", "\n", " # change alpha values\n", " color_array[:,-1] = np.linspace(0.0, 1.0, ncolors)\n", "\n", " # create a colormap object\n", " map_object = LinearSegmentedColormap.from_list(name=newmap, colors=color_array)\n", "\n", " # register this new colormap with matplotlib\n", " plt.register_cmap(cmap=map_object)" ] }, { "cell_type": "code", "execution_count": 3, "id": "65730d47-362c-4378-851a-2a2bdb295a09", "metadata": {}, "outputs": [], "source": [ "mymesh = xr.open_dataset('/home/sallen/MEOPAR/grid/mesh_mask201702.nc')\n", "depthf = interp.interp1d(mymesh.z, mymesh.gdepw_1d)\n", "figrid = 386\n", "igrid = figrid-1\n", "fjgrid = 177\n", "jgrid = fjgrid-1\n", "tmask = mymesh.tmask[0]\n", "fmask = mymesh.fmask[0, :, igrid]\n", "fmask2 = mymesh.fmask[0, :, :, jgrid]\n", "f_lons = mymesh.glamf[0, igrid]\n", "f_lats = mymesh.gphif[0, :, jgrid]\n", "u_lons = mymesh.glamv[0, igrid] # note switch to v with the switch from t to f\n", "v_lats = mymesh.gphiv[0, :, jgrid] \n", "w_depths = mymesh.gdepw_1d\n", "xs, ys = np.meshgrid(np.array(mymesh.glamv[0, igrid]), np.array(mymesh.gdept_1d))" ] }, { "cell_type": "code", "execution_count": 4, "id": "bd195fde-b909-4e5e-8193-712953cf5fcf", "metadata": {}, "outputs": [], "source": [ "m2lon = ((u_lons[314]-u_lons[260])/(mymesh.e2f[0, igrid, 260:314].sum())).values\n", "m2lat = ((v_lats[301]-v_lats[234])/(mymesh.e2f[0, 234:301, jgrid].sum())).values" ] }, { "cell_type": "code", "execution_count": 5, "id": "ba51d77b-c01b-4e91-b004-5095b386546a", "metadata": {}, "outputs": [], "source": [ "def draw_patches(fmask, ll, ul, u_lons, w_depths):\n", " maxw = w_depths.shape[0] - 1\n", " topo = []\n", " for i in range(ll,ul):\n", " for j in range(maxw):\n", " if fmask[j, i] == 0:\n", " rect = Rectangle((u_lons[i], w_depths[j]), u_lons[i+1]-u_lons[i], w_depths[j+1]-w_depths[j])\n", " topo.append(rect)\n", " pc = PatchCollection(topo, facecolors='lightgray', edgecolors='lightgray')\n", " return(pc)" ] }, { "cell_type": "code", "execution_count": 6, "id": "f2adb6a9-8f22-48e6-9b0f-5c9c687a4e89", "metadata": {}, "outputs": [], "source": [ "def draw_patches_vic(fmask, ll, ul, v_lats, w_depths, alpha=1):\n", " topo = []\n", " for i in range(ll, ul):\n", " for j in range(39):\n", " if fmask[j, i] == 0:\n", " rect = Rectangle((v_lats[i], w_depths[j]), v_lats[i+1]-v_lats[i], w_depths[j+1]-w_depths[j])\n", " topo.append(rect)\n", " pc = PatchCollection(topo, facecolors='lightgray', edgecolors='lightgray', alpha=alpha)\n", " return(pc)" ] }, { "cell_type": "code", "execution_count": 7, "id": "d3512bee-77bb-4980-8162-155ea1e192c8", "metadata": {}, "outputs": [], "source": [ "def get_data_forward(month, year, section=2):\n", " amonth = glob.glob('/data/sallen/results/Ariane/FullNorth/*'+month+str(year)+'/ariane_positions_quantitative.nc')\n", " mydata = xr.open_dataset(amonth[0])\n", " lons = mydata.init_lon[(mydata.final_section==section)]\n", " depths = depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]\n", " transports = mydata.init_transp[(mydata.final_section==section)]\n", " sals = mydata.init_salt[(mydata.final_section==section)]\n", " mydata.close()\n", " for f in amonth[1:]:\n", " mydata = xr.open_dataset(f)\n", " lons = np.concatenate((lons, mydata.init_lon[(mydata.final_section==section)]))\n", " depths = np.concatenate((depths, depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]))\n", " transports = np.concatenate((transports, mydata.init_transp[(mydata.final_section==section)]))\n", " sals = np.concatenate((sals, mydata.init_salt[(mydata.final_section==section)]))\n", " mydata.close()\n", " return lons, depths, transports, sals" ] }, { "cell_type": "code", "execution_count": 8, "id": "544f4a8a-62e8-48e1-b0cf-8afe18064633", "metadata": {}, "outputs": [], "source": [ "def get_data_forward_vic(month, year, section=2):\n", " amonth = glob.glob('/data/sallen/results/Ariane/FullSouth/*'+month+str(year)+'/ariane_positions_quantitative.nc')\n", " mydata = xr.open_dataset(amonth[0])\n", " lats = mydata.init_lat[(mydata.final_section==section)]\n", " depths = depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]\n", " transports = mydata.init_transp[(mydata.final_section==section)]\n", " sals = mydata.final_salt[(mydata.final_section==section)]\n", " mydata.close()\n", " for f in amonth[1:]:\n", " mydata = xr.open_dataset(f)\n", " lats = np.concatenate((lats, mydata.init_lat[(mydata.final_section==section)]))\n", " depths = np.concatenate((depths, depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]))\n", " transports = np.concatenate((transports, mydata.init_transp[(mydata.final_section==section)]))\n", " sals = np.concatenate((sals, mydata.final_salt[(mydata.final_section==section)]))\n", " mydata.close()\n", " return lats, depths, transports, sals" ] }, { "cell_type": "code", "execution_count": 9, "id": "d26df1c3-cb73-4964-9ce9-095441838b80", "metadata": {}, "outputs": [], "source": [ "def get_data_backward_gi(month, year, section = 2):\n", " amonth = glob.glob('/data/sallen/results/Ariane/InGIslands/*'+month+str(year)+'/ariane_positions_quantitative.nc')\n", " mydata = xr.open_dataset(amonth[0])\n", " lons = mydata.init_lon[(mydata.final_section==section)]\n", " depths = depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]\n", " transports = mydata.init_transp[(mydata.final_section==section)]\n", " sals = mydata.init_salt[(mydata.final_section==section)]\n", " mydata.close()\n", " for f in amonth[1:]:\n", " mydata = xr.open_dataset(f)\n", " lons = np.concatenate((lons, mydata.init_lon[(mydata.final_section==section)]))\n", " depths = np.concatenate((depths, depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]))\n", " transports = np.concatenate((transports, mydata.init_transp[(mydata.final_section==section)]))\n", " sals = np.concatenate((sals, mydata.init_salt[(mydata.final_section==section)]))\n", " mydata.close()\n", " return lons, depths, transports, sals" ] }, { "cell_type": "code", "execution_count": 10, "id": "3a06a06a-dabc-491c-9ef6-fe6899a27eb0", "metadata": {}, "outputs": [], "source": [ "def get_data_backward_gi_vic(month, year, section = 3):\n", " amonth = glob.glob('/data/sallen/results/Ariane/BackSouth/*'+month+str(year)+'/ariane_positions_quantitative.nc')\n", " mydata = xr.open_dataset(amonth[0])\n", " lats = mydata.init_lat[(mydata.final_section==section)]\n", " depths = depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]\n", " transports = mydata.init_transp[(mydata.final_section==section)]\n", " sals = mydata.final_salt[(mydata.final_section==section)]\n", " mydata.close()\n", " for f in amonth[1:]:\n", " mydata = xr.open_dataset(f)\n", " lats = np.concatenate((lats, mydata.init_lat[(mydata.final_section==section)]))\n", " depths = np.concatenate((depths, depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]))\n", " transports = np.concatenate((transports, mydata.init_transp[(mydata.final_section==section)]))\n", " sals = np.concatenate((sals, mydata.final_salt[(mydata.final_section==section)]))\n", " mydata.close()\n", " return lats, depths, transports, sals" ] }, { "cell_type": "code", "execution_count": 11, "id": "e0128b13-d5af-4c5d-8478-2bd31897f43e", "metadata": {}, "outputs": [], "source": [ "def get_data_forward_gi(month, year, section = 2):\n", " amonth = glob.glob('/data/sallen/results/Ariane/SouthGIslands//*'+month+str(year)+'/ariane_positions_quantitative.nc')\n", " mydata = xr.open_dataset(amonth[0])\n", " lons = mydata.init_lon[(mydata.final_section==section)]\n", " depths = depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]\n", " transports = mydata.init_transp[(mydata.final_section==section)]\n", " sals = mydata.init_salt[(mydata.final_section==section)]\n", " mydata.close()\n", " for f in amonth[1:]:\n", " mydata = xr.open_dataset(f)\n", " lons = np.concatenate((lons, mydata.init_lon[(mydata.final_section==section)]))\n", " depths = np.concatenate((depths, depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]))\n", " transports = np.concatenate((transports, mydata.init_transp[(mydata.final_section==section)]))\n", " sals = np.concatenate((sals, mydata.init_salt[(mydata.final_section==section)]))\n", " mydata.close()\n", " return lons, depths, transports, sals" ] }, { "cell_type": "code", "execution_count": 12, "id": "c1e92765-664e-44c4-a453-c19c6fb78db9", "metadata": {}, "outputs": [], "source": [ "def get_data_forward_gi_vic(month, year, section = 3):\n", " amonth = glob.glob('/data/sallen/results/Ariane/FullSouth//*'+month+str(year)+'/ariane_positions_quantitative.nc')\n", " mydata = xr.open_dataset(amonth[0])\n", " lats = mydata.init_lat[(mydata.final_section==section)]\n", " depths = depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]\n", " transports = mydata.init_transp[(mydata.final_section==section)]\n", " sals = mydata.final_salt[(mydata.final_section==section)]\n", " mydata.close()\n", " for f in amonth[1:]:\n", " mydata = xr.open_dataset(f)\n", " lats = np.concatenate((lats, mydata.init_lat[(mydata.final_section==section)]))\n", " depths = np.concatenate((depths, depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]))\n", " transports = np.concatenate((transports, mydata.init_transp[(mydata.final_section==section)]))\n", " sals = np.concatenate((sals, mydata.final_salt[(mydata.final_section==section)]))\n", " mydata.close()\n", " return lats, depths, transports, sals" ] }, { "cell_type": "code", "execution_count": 13, "id": "c0ad7b16-d9c0-4ea9-aced-d733fb9855c4", "metadata": {}, "outputs": [], "source": [ "def get_data_backward(month, year, section=2):\n", " amonth = glob.glob('/data/sallen/results/Ariane/BackNorth/*'+month+str(year)+'/ariane_positions_quantitative.nc')\n", " mydata = xr.open_dataset(amonth[0])\n", " lons = mydata.init_lon[(mydata.final_section==section)]\n", " depths = depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]\n", " transports = mydata.init_transp[(mydata.final_section==section)]\n", " sals = mydata.init_salt[(mydata.final_section==section)]\n", " mydata.close()\n", " for f in amonth[1:]:\n", " mydata = xr.open_dataset(f)\n", " lons = np.concatenate((lons, mydata.init_lon[(mydata.final_section==section)]))\n", " depths = np.concatenate((depths, depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]))\n", " transports = np.concatenate((transports, mydata.init_transp[(mydata.final_section==section)]))\n", " sals = np.concatenate((sals, mydata.init_salt[(mydata.final_section==section)]))\n", " mydata.close()\n", " return lons, depths, transports, sals" ] }, { "cell_type": "code", "execution_count": 14, "id": "6f9f9a42-b3de-4bac-8157-cedce7f662ee", "metadata": {}, "outputs": [], "source": [ "def get_data_backward_vic(month, year, section=2):\n", " amonth = glob.glob('/data/sallen/results/Ariane/BackSouth/*'+month+str(year)+'/ariane_positions_quantitative.nc')\n", " mydata = xr.open_dataset(amonth[0])\n", " lats = mydata.init_lat[(mydata.final_section==section)]\n", " depths = depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]\n", " transports = mydata.init_transp[(mydata.final_section==section)]\n", " sals = mydata.final_salt[(mydata.final_section==section)]\n", " mydata.close()\n", " for f in amonth[1:]:\n", " mydata = xr.open_dataset(f)\n", " lats = np.concatenate((lats, mydata.init_lat[(mydata.final_section==section)]))\n", " depths = np.concatenate((depths, depthf(mydata.init_z[(mydata.final_section==section)]-1.)[0]))\n", " transports = np.concatenate((transports, mydata.init_transp[(mydata.final_section==section)]))\n", " sals = np.concatenate((sals, mydata.final_salt[(mydata.final_section==section)]))\n", " mydata.close()\n", " return lats, depths, transports, sals" ] }, { "cell_type": "code", "execution_count": 15, "id": "f4eea84b-7986-45fb-bab2-869dec8e9df9", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "15\n", "16\n", "17\n", "18\n" ] } ], "source": [ "month = {1:'jan', 2: 'feb', 3: 'mar', 4: 'apr', 5: 'may', 6: 'jun', 7: 'jul', \n", " 8: 'aug', 9: 'sep', 10: 'oct', 11: 'nov', 12: 'dec'}\n", "lons, depths, transports, sals = np.array([]), np.array([]), np.array([]), np.array([])\n", "months = 4*np.ones_like(lons)\n", "years = 15*np.ones_like(lons)\n", "for y in [15, 16, 17, 18]:\n", " print(y)\n", " for m in month:\n", " newlons, newdepths, newtransports, newsals = get_data_backward(month[m], y)\n", " newmonths = m*np.ones_like(newlons)\n", " newyears = y*np.ones_like(newlons)\n", " lons = np.concatenate((lons, newlons))\n", " depths = np.concatenate((depths, newdepths))\n", " transports = np.concatenate((transports, newtransports))\n", " sals = np.concatenate((sals, newsals))\n", " months = np.concatenate((months, newmonths))\n", " years = np.concatenate((years, newyears))\n", " \n", " newlons, newdepths, newtransports, newsals = get_data_backward_gi(month[m], y)\n", " newmonths = m*np.ones_like(newlons)\n", " newyears = y*np.ones_like(newlons)\n", " lons = np.concatenate((lons, newlons))\n", " depths = np.concatenate((depths, newdepths))\n", " transports = np.concatenate((transports, newtransports))\n", " sals = np.concatenate((sals, newsals))\n", " months = np.concatenate((months, newmonths))\n", " years = np.concatenate((years, newyears))" ] }, { "cell_type": "code", "execution_count": 16, "id": "58993599-6295-4624-827f-9d5c9acedc9b", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "15\n", "16\n", "17\n", "18\n" ] } ], "source": [ "month = {1:'jan', 2: 'feb', 3: 'mar', 4: 'apr', 5: 'may', 6: 'jun', 7: 'jul', \n", " 8: 'aug', 9: 'sep', 10: 'oct', 11: 'nov', 12: 'dec'}\n", "olons, odepths, otransports, osals = np.array([]), np.array([]), np.array([]), np.array([])\n", "omonths = 4*np.ones_like(olons)\n", "oyears = 15*np.ones_like(olons)\n", "for y in [15, 16, 17, 18]:\n", " print(y)\n", " for m in month:\n", " newlons, newdepths, newtransports, newsals = get_data_forward(month[m], y)\n", " newmonths = m*np.ones_like(newlons)\n", " newyears = y*np.ones_like(newlons)\n", " olons = np.concatenate((olons, newlons))\n", " odepths = np.concatenate((odepths, newdepths))\n", " otransports = np.concatenate((otransports, newtransports))\n", " osals = np.concatenate((osals, newsals))\n", " omonths = np.concatenate((omonths, newmonths))\n", " oyears = np.concatenate((oyears, newyears))\n", " \n", " newlons, newdepths, newtransports, newsals = get_data_forward_gi(month[m], y)\n", " newmonths = m*np.ones_like(newlons)\n", " newyears = y*np.ones_like(newlons)\n", " olons = np.concatenate((olons, newlons))\n", " odepths = np.concatenate((odepths, newdepths))\n", " otransports = np.concatenate((otransports, newtransports))\n", " osals = np.concatenate((osals, newsals))\n", " omonths = np.concatenate((omonths, newmonths))\n", " oyears = np.concatenate((oyears, newyears))" ] }, { "cell_type": "code", "execution_count": 20, "id": "5c927e83-dca2-4fe0-9645-1a85891f8a6b", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "oo 31870.904511511468 26740.307183473386 5130.597328038051\n", "11.059173899088456\n", "oo 25952.223241768603 19524.70486665357 6427.518375115036\n", "2.643509514932847\n", "3.4226393475352457\n", "0.761932724039425\n", "45212.61098342584\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA6oAAAJtCAYAAAA/5o2mAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABNY0lEQVR4nO39ebhlZ1kn/H/vGlKVylgZyTxABhIig0EEIQQbERXhbRG7AVFQAW3ttqVfBxyAX9v6anfT3bZ0tyJ0o6i8gDjw0oJKS5gESZhDyEhC5rlSqarUXM/vj30qWbVzhn32PufUOiefz3XtK2s967mf/ZzKOXvv715TtdYCAAAAfbHqYE8AAAAAugRVAAAAekVQBQAAoFcEVQAAAHpFUAUAAKBXBFUAAAB6RVAFAACgVwRVAAAAeuUxF1Sr6vSq+pWq+mxV3V5VO6vqlqr6eFX9bFUdf7DnCAAA8FhWrbWDPYclU1X/MslvJzl0lm73J/mp1tr7lmZWAAAAdD1mgmpV/XqSXx1qvi7J7UlOTfL4oW0/1lr7X0sxNwAAAB7xmAiqVfUDST7Qaboqyataa1/o9Lk4yR8leeJU054k39Fa+9ySTRQAAICVH1Sram2Sq5OcPdV0a5Jvaa1tmqbvMUm+kuSUqaZPttYuWZKJAgAAkOSxcTGlV+SRkJokb5gupCZJa+3+JG/oND2nqgRVAACAJfRYCKov6yzfnuQv5uj/51P9pqsHAABgka3ooFpVhyZ5fqfpI621PbPVTG3/m07TixdjbgAAAExvRQfVJBckWddZ//SIdd1+p0+duwoAAMASeCwE1a7rRqwb7jc8DgAAAItkpQfVM4fWbx6x7ptzjAMAAMAiWXOwJ7DIjhxaf2DEus1D60fM1LGqXpfkdUly6KGHfutZZ5018uT6YF/bl6tvveHh9dWrJvvu4rANG3L2Sacf0LZrx+6JxuSxo6qyfv36TktLUmOPd8P9t0xQnTxw3/BLwfwdsmZN9k1wG7C1q9dk997J/obWrlmbPXtnPT1/Vocdsi67J6g/7rD1OWLDjC+j09qb9XN3YkHsa/ty3R03jdT3xHVr09reiZ5vw9q12bx7/DHWrlmbTTt3TTSH1atWZfuu8ceoquzbt2+i+klvDzjxGJVMeofCqgxepseun+xnWLVq1cT/jqtqVfa18f9frlm1Knv3TTaHNatXZ+++Cf4mVq/Jngnqk8F71Z6944+xbs3a7JrgfeLQtWuzc89k73Ub1h6SHbvHH+PwQ9Zm+wT1SXL02lXZvmf8f4c1D26f95/U0Weel5rw8/tCu+qqq+5trR0/6TgrPagePrS+fcS64X4zfsJqrb09yduT5MILL2zvfe97R59dT1z0cy98ePmIw4b/yebnWd/y5PzJr/z3A9puve6OicbkseXCCy98eHlf25uq8aPmD/3JGyYKqh96z/+eoHrgtGNPzI7d438gPuGoY3L7/XdPOIcTctfm+8eu/7Yzz8pdm+4bu/4nnnFennfxpfOqecAZF0vq0je/fKR+P/2E47Nj+9aJnutppz4uf33bA2PXn3LcSfmzb9wy0RwOW39orrrt1rHrD1mzNlt2bBu7fu2aNdk14QfitWvWZPckH4hXr56ofv8cJgo3aw+ZaA4b1q3PzgleX5PkiPUb8tCuHWPXbzzs8GzZPurHy+kdd8RReeCh8f+uTj76mNy39cGJ5nDGscfnzs3T3r1xJOec+Ljcct+9Y9dfePIpufGeu8auT5KnnnJqrr7rzrHrn3Xq4/L1CeqT5EUnHJYr7xn/3+GYj301bZ5/Uy993/vHfr7FctFFFw0fnTqWfsXvhbd2aH3UV8Phd4/hcZaNub5pnPSbSIDF4bVpqXgfAFi+VvJr+EoPqsNfd456LNmhc4yzbMy1N2qSvVUAi8dr01LxPgCwfK3k1/CVHlSHj6PYMGLdcL8tCzCXg8IeVWB58tq0VLwPACxfK/k1fKUH1XuG1k8asW643/gHmx9k9qgCy5PXpqXifQBg+VrJr+ErPahePbR+xoh1w/2Gx1k27FEFlievTUvF+wDA8rWSX8NXelD92tD600as6/bbleT6hZnO0rNHFVievDYtFe8DAMvXSn4NX9FBtbV2S5JvdJqeO2Jpt9+n2qQ3jTuI7FEFlievTUvF+wDA8rWSX8NXdFCd8hed5Uur6vTZOk9t7wbVDyzKrJaIParA8uS1aal4HwBYvlbya/hjIaj+ryT7ppZXJfm1Ofq/KY/8u2xN8r5FmhcAAADTWPFBtbX2tSR/3Gn6iar6ien6VtXrk/x4p+k/ttaW7RV/E4f+AsuV16al4n0AYPlaya/haw72BJbIzyd5TpKzptb/oKq+P8n/m+T2JKckeXmSF3VqLk/yH5ZykovBob/A8uS1aal4HwBYvlbya/hjIqi21u6uqu9N8jdJ9p+j+uKpx3S+kuRFrbWHlmJ+i6m1Nusv8Er+FgZYzlqE1aXhfQBg+Zrrs/5ytuIP/d2vtXZ1kouS/PcMzj2dzn1J/l2Sp7fW7l6quS0me1SB5clr01LxPgCwfK3k1/DHxB7V/VprDyb56ar6v5NcmuSMJBuT3JvkhiSfaK3tOXgzXHj2qLLcfO1rw7c/fsSqVaty/hPPW8LZcPDYo7pUvA8ALF8reY/qYyqo7tda257kwwd7HkvBHlVWkn379s3diRXCa9NS8T4AsHyt5Nfwx2RQfSyxR5WV5qqvfX3Gbfa4riT2qC4V7wMAy9dK3qP6mDlH9bHKHlUeS+xxXUm8Ni0V7wMAy9dKfg0XVFc491EFlievTUvF+wDA8rWSX8MF1RXOHlVgefLatFS8DwAsXyv5NVxQXeHsUQWWJ69NS8X7AMDytZJfwwXVFc4eVWB58tq0VLwPACxfK/k1XFAFAACgV9yeZoVzexpgeXJ7mqXSWstl/7/3zLh9287t+b7f/LElnBEAo3J7GpYth/4Cy5PXpqUy1/vAYesOXaKZADBfK/mzvD2qK5w9qjzWXPW1rz+8/JanvPZR20f5m9i//c1j1G/b8VAuffMr5jNlpmWP6lLxPgGwfNmjyrJljyocaNK/iTn3Pq3fMO85MR2vTUvF+wTA8rWSX6MF1RXO7WngQJP+TfibWir+HZeK32mA5Wslv0YLqiucb8rhQIu9R9Xf1ELx77hU/E4DLF8r+TVaUF3hfFMOB7JHdbnw77hU/E4DLF8r+TVaUF3hfFMOB7JHdbnw77hU/E4DLF8r+TXaVX/JiUdtfHj51eefli3bt4891tNOPS7bdzywALOC5esJp5728PL5T3hh9uzbM/ZY55x8Wr5y5Ucmms/L77g9D955x9j1u26/Pzu23D92/S1XfCnv3PbOedW89H2fH/v5WHinHHtikuSfrNqTfav3TTTWabU2Rx69duz6o086KVftmeyWOaccd3K27B6//vTjTsl1d9w4fv3xJ+fmu28dfwJJTjvhpNxyz/h/16cdf1Juufu2ieZw9omn5eZ7bh+7/rxTzsg37xq/PknWr53so+y3Pu2CXH3P+P8vXnD+k3L1hP8vLz39zFx5xy1j1z/35FPy5Qn/HZ97zNH5wr33jF1/yeGH5UubN87dcQbP2nBovrTlkLHrk+TiQ9bmquOOGbv+ybUn162f7LXl9NvuypN2bBm7/hNHH5edm8d/v11pBNUVbr63Hdi+a1ce2rlj7OfbsmN7Dl17+FDr+MEXlqPL/9uHH17+3l9/y0Rj3bt5U+7b8sBEY+zcujmrJwjLO7Zsy54dk/0dt33zCzcr+XL7y9Hn/8fgy5L3vvaitEwWVA8/8sFc/8CDY9cfs+Ge3HD3fRPNYduelm0T/E5v3r41u/eMn3R37NyR3Xv3jl2fJDt37snuPeP/Xe/atTt7JpzDjl07Jxpjx86d2bV3/H/H8046Ld+4a7KwvXn79hx+2BFj12/bszOPP+XMieawu/blvDPPHbt+755tOf/M8yaaw65td+UJZ5wzdv3OTbfltJPPGLv+oU235/jjHzd2/WCMO3PkMSeMX3/f7TnsmOMnmsPOW+/KIRuPG7t+3+5veL/scOjvCueQLmA58toEAHNbye+XguoK5yIZwHLktQkA5raS3y8F1RXOHlVgOfLaBABzW8nvl4LqCmePKrAceW0CgLmt5PdLQXWFs0cVWI68NgHA3Fby+6WgCgAAQK8IqiucQ3+B5chrEwDMbSW/XwqqK5xDf4HlyGsTAMxtJb9fCqoAAAD0iqAKAABAr6w52BNgcbXWZj0kYCUf1w4Hy/VXX//w8n992Q8fsG3H7t35hb9871JPadmZ67WLpXX1169Okjz5De9/1La9Ox/KV//bq5Z6SgBkZb9f2qO6wjlHFfpl/dq1B3sKy4LXpuVj9boNB3sKAI9ZK/n9UlAFAACgVwRVAAAAekVQBQAAoFdcTGmFczEl6Bd/cwMv+aNPZe2Gw2bcvpIvDrHS+J0GOHhW8vulPaornIspQb/4mxuYLaQm/p2WE/+vAA6elfwaLKgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArguoKN9dtA9xWAJaWv7kBr00rh/9XAAfPSn4NFlRXOLengX7xNzfgtWnl8P8K4OBZya/BgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvbLmYE+AxdVaS1XNuv053/L4h9dfuWFH7tjcxn6+0w7dnT27t4xdDytday337nzayP1r8+E55fivHND2a//XT2fd2kNmfY7Z/u737tqZy3//rSPP4dKzT8yeqz8ycv9hJ5x2Sr5+w53zmuNc2+mP1lpOe/oLRu6/5+xn5oR7rh/7+dae+JRccMxX5u44i3OOeWIOOfzysesvPPviXH3TF8euv+icb83V3/zS2PVJ8rQLvy1X3/Hlseu//QnfmuvuuHKiOTzrrG/J1TdPMIfTL8yXbxi/PkleuOvcieqfedqpueb+W8euf9qGw3PXnm0TzeHs7dtz77YHx64/YfOmPDjhrqfDH9icrROMsX7brmxvN4xdv2bvuhy6b/z/D0mSfYdn446vj11+Tx2XnZuum2gKX919UjZ/c/wxnvfE47Pn3vn9Pq3k90tBdYWb6xe3qrK30+eaB7bk1q3bx36+XVu357y1Rw217hh7PFhpqiqnH3/syP2PWL827/7Jtx/QdvXXZ38jnuvvfvUh67J3166R59B278izfuznRu4/bN+OB3LOj19yQNvVW+Z+bWJ5qKpsfd73j9z/izffmzs2DL9PzMPme7Np19bx65PcdPetuX/rprHr79l8T+5/cPz6e7fcny07JvsZ7tt2fw5Zv27s+q17tuXcM86eaA5Z2/Lcpzxz7PL1e/bmJc9+4dj1q7duyoYNR45dnyS5/5acfco5Y5e3O6/JyUc8YaIp7Lv5yznq6JPGrt+z9YoctvGUieawY9OXc+xRp41d/9Cmq3P4xvF/hp133ZoNG08Yuz5Jtt15R9YdfvTY9bvvvDurVs38JfAoHrrn9uye4EuHvVsfyjnf8V3zqrl7Bb9fOvQXAACAXhFUAQAA6BVBFQAAgF4RVAEAAOgVQRUAAIBeEVQBAADoFUEVAACAXhFUAQAA6BVBFQAAgF4RVAEAAOgVQRUAAIBeEVQBAADoFUEVAACAXhFUAQAA6BVBFQAAgF4RVAEAAOiVNQd7Aiyu1lqqatbtwNLpw9/cXHN4xs+8MWvWrTug7ZpHD5LM8tryqO1bHtUhySz1c26nL1oO/u80wEp1zwt/OW3Nupk7zPV+vIzZo7rCzRZSR9kOLKw+/M3NNYfhkDrDIJNtnzOEHvx/J0ZT/l8BLJpZQ2qyYkNqIqgCAADQM4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL3Sm6BaVWdWVRvjcfU8n+f8qvrNqvp8Vd1dVTuq6qaq+puq+omqOmKxfkYAAADmtuZgT2CpVNWaJG9O8sYkq4c2nzH1eEGSX6uqV7fWPrbEUwQAACD9DqqfSLJ9hH63jDjeO5P8SGe9Jfl6knuTnJXktKn205P8bVV9X2vtb0ccGwAAgAXS56D6o621mxZioKp6Qw4MqZ9I8trW2rWdPs9P8odJTs7g3+X9VfUtrbVvLsQcAAAAGE1vzlFdLFV1bJI3dZq+mOQF3ZCaJK21jya5JMnWqaYjk/z6kkwSAACAh634oJrkZ5Ic1Vl/fWtt53QdW2s35MBw+sqqOnMR5wYAAMCQx0JQfVln+XOttcvn6P+OJDumllcleemizAoAAIBpreigWlVnJ7mw0/ShuWpaa/cn+Uyn6cULPS8AAABmtqKDapKnDq1/esS6br+nLMxUAAAAGEWfg+q/r6qvVtUDVbWzqu6sqsur6r9U1XNHHOOCofXrRqzr9juyqk4dsQ4AAIAJ9TmovizJkzK4ENIhSU5McnGSn01yWVX9Y1U9aY4xzuws701y+4jPPXxLmjOn6wQAAMDC63NQvS/J55L8nyT/mOTeoe3fluRzVfWiWcY4srO8pbW2d8Tn3jy0fsRMHavqdVV1RVVdsWnTphGHBwAAYCZ9C6qfT/LTSc5urR3XWntGa+35rbVvb60dn8Ee1T/r9D80yXur6ikzjHd4Z3n7POYx3HfGoNpae3tr7eLW2sUbN26cx1MAAAAwnTUHewL7tdZuyiCIztbn80leVlX/Msl/nWreMLV8yTQlazvLe+YxneG+a6fttUJ82xPPenj527cfmVsfuGfssR534kmpnVsWYlqwYm3cOON3X4+ybs3ifJ/4+IuGT+GH8a2rY0buu7OOy+pDjh/7udasPi7nn3bc2PUDp2TtkaePXX3ohjPyxCNOGbt+/WGn54mHnzx2fZIctuHUHHPY3WPXH77qhGzcN9mRYOu2H561u+8fu752r86qe68fu74l2XPN58euT5Kd996fXTuGD9ob3db7t2b1zgcnmsN92/Ym28cfY9Oe9WnbJ/t32H7ocdn74KfGrl9z3KnZc83Xxq6vxz0hO2+6aez6JFl36nl56Prxf5/2nXhONt/49YnmcNSZ56XtHfUAzke7ZsPjcu7hR000h5WkN0F1Plprvzu1F/XHppqeU1UXt9auGOq6rbO8fh5PMdx327S9loHWWqpq1u2nnHXRw+tf+MYV2Xfq+B8g9m04OuetO3qo9aGxx4OVpiX5xW/bMGef7l/t1V+f7I1zOo97wQ8t+Jg8dj1u9yxvsVU59tRHQuF199+c444+aeznumfb3nzrmZN90fLV27flgpPPG7v+xvt25wmnnzN2/c33784TTjt77PokufvB3Tn/2PHD9uatD+Wso46daA4PPfhgjj/08Lk7zmDn9rtzxPpDxq7fftv1Wb92so+yex68K6vXrRt/gPvvSB129ERz2HvfrVl12PjhZOf9t2f1+sMmmsOOu27OqjXj75fZeduNWT3BF6u7br0+bc+useuTZNc3r8m+HeN/ZN+6/epsv+/OiebQ9u3NjgfH/wJo/aHrcszrP3hA253XjB++l7u+Hfo7H78xtP7Cafps7SzP/snwQMN9l+0uwtlC6ijbgYU1yl+cv0pWlNYO9gwAWIaWbVBtrX0jB16d9/xpunWPYT2sqkY93m74q97xjwkBAABgXpZtUJ1yR2d5uhNWrh5aP2PEcbv99iW5dj6TAgAAYHzLPah2D9Gd7qq+w2d1P23Ecbv9bmqtzeeKwQAAAExg2QbVqlqX5AmdpunOfr4iB14I6bkjDt+9gvBl85sZAAAAk1i2QTXJD+TAPaqPuqb21J7Qj3SaXlpVs15UqaqenaR7Kb4PTDJJAAAA5mdZBtWqOiHJb3WatuXAQNr1zs7yUUl+bo7h39xZvjnJR+c9QQAAAMbWi6BaVc+sqt+rqjlvalZVFyX5WJLuzcPe2lq7b7r+rbUPJ/l4p+lNVfW9M4z9G0me3+3bWpvspk4AAADMy2R3SV4465K8Psnrq+rLSf4+yVcyOO90S5LDMzgf9buTfF8ODNh/l0ffU3XY65J8JskxSQ5J8sGqek+Sv0xyX5KzkrwmyXM6NR9M8u5JfigAAADmry9BtevJU49R/GGSn55rr2dr7dqqekmSv8ogrK5O8sNTj+n8fZKXt9b2jTgPAAAAFkgvDv1NcmOSP05y/Qh992Swt/OftNZe3VrbNldBkrTWPpXkwiTvSbJzhm63JnlDku9qrT00yrgAAAAsrF7sUW2tfTPJq5Kkqo5L8i1Jjk9ybJKNSXYk2ZTkuiRXjHtf09banUleUVVHJbk0yalJjkhyV5JrknymtdYm+mEAAACYSC+Caldr7d4MDr1dzOfYnMFhwAAAAPRMXw79BQAAgCSCKgAAAD0jqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL2y5mBPgIPv6LWrH14+96hjs2nXtrHHOmHd2rRdDy7EtABYIXbuag8vX3blluxr+8Yeq9asza7dk318uf2Byp49489h07bVuW3VnrHrd+2qbL1/89j1SbJ63yF56JZbxq7ft/qoPHDDlRPNYdvhp2bL9f8wdv32o8/OA1/46ERz2H7/XRPVH3HSGbn/6k+PXb/hcadny1evnmgOhx17Urbc9M2x6w8/8bQ8dOftE83hyFMfn6133TZ2/cYzzsmDt900dv1x539LHrzlG2PXJ8lxT/zW3H/D18auP+XcJ+fea7440RyOOev83HP1+GO0ubs8pgiqK15LUrNuP/+Upz68tm3Pjhx5yOHjP12tTq07cqhx6/jjwQrTWlKz/UmO2AeWk+MOP/Th5ZOOODRrVo3/C/7gnpbHbzxqsgnteygXnnjY2OU33LslTzxu3dj199x3by48YvygmyT333N7zsymsevvu/PrOW7rZOHmpm98LXXv+AHr/quvytabrxm7fs0Rx2bbnTeNXZ8kuzbdle333T12/Z4dD2Xn1sm+dNi3a1d2bRl/jIdWrcqubVsmmsO2u27Nrs0PjF2/9fZvZsfm+8au33LLDRP9fxiMcX12TPDFxZZbb8iO+yecwyHrs2PzvSP3v/Q3/jBr1m84oO2qa66faA4riUN/V7y5Pgz4NAxLaZQAKqQCwMo3HFI5kKAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPTKmoM9AQ6+1Z3lY7I62bN77LE2rKq0XdsmnxQAK9KXb96cqvHrt+7cmwe27JxsEmtXZVXbM3b5rl27c/9tt45dv23vYfnc5z45dn2SbDz3Kbn5C+OPcdz5F+e6b14z0RyOv/AZ2bRp/Pf8k7/j2bn3qM9NNIcNpz95ovqjzn1G7rvqH8eu33j+07PpmvHrk+TkJ16cHdd+euz60y781uy9/h8mmsPRT/r2bLv+8rHrj3nSM3Pv1eP/vzz+om/Pvdd+aez6JNlwxtPz0PVfGLv+npMvTY778kRzOOOcJ2T3l/5iojF4hKBKztj4xIeXd+5bnarxd7TX3p1Zte6oodbNY48HwPJ30w2PBKLf/d5jDti2a2/LO786emi85q5tOe3o9RPN56F9e3PusYeMXb/lzjty4aq7x67/x6/cmL1fuWzs+iTZfc8N2XHN+MFi233fzLZbvjHRHNbv2ZWtt48/xrojjs/2++4cu/7Ic749275x5dj1SfLQpgeye/v2set3PLgtux8avz5J9mzbmrWHrBu7ft/2bTnulNMnmkPbsyOnPvGisevrkORx3/Wiset3bn0wxz/n+WPXJ8ntX7s+G48e/gw6umtv/Wa23jP+72OSbNuwOkesHv2buNZaapJv7lY4h/4CAAfNIfP4UAewkgipsxNUAQAA6BVBFQAAgF4RVAEAAOgVQRUAAIBeEVQBAADoFUEVAACAXhFUAQAA6BVBFQAAgF4RVAEAAOgVQRUAAIBeEVQBAADoFUEVAACAXhFUAQAA6BVBFQAAgF4RVAEAAOiVeQfVqjqyqi6tqn9TVe+pqmural9VtanHZZNOqqrOr6rfrKrPV9XdVbWjqm6qqr+pqp+oqiMmGPviqvqdqvpKVd1XVQ9V1fVV9ZdV9fKqWjfp/AEAABjfmvl0rqprkpyTpBZjMlW1Jsmbk7wxyeqhzWdMPV6Q5Neq6tWttY/NY+zDkvynJK+bZvPjpx4vSXJlVb2qtfal+f8EAAAATGq+e1TPzSKF1CnvTPKreSSktiRXJflEkls6/U5P8rdV9YJRBq2qtUk+mAND6u4kX0nyySR3ddqflOQTVfXkcX4AAAAAJjPuOapbMgiP/znJDyf54qQTqao3JPmRTtMnkpzfWruwtfbc1trpSb4rye1T29ckeX9VnTHC8G9N8p2d9Q8kOau19uTW2iVJTk7yz5I8OLX9iCQfmuQQYwAAAMYz36D6yiTnJzlqKjy+obX2J3kk4I2lqo5N8qZO0xeTvKC1dm23X2vto0kuSbJ1qunIJL8+x9jnJfmpTtOHkrystXZbZ9x9rbX3JfmeJHunmk9N8vPz/2kAAACYxLyCamvtT1tr17TW2gLP42eSHNVZf31rbecMc7ghB4bTV1bVmbOM/Ut55Fzc3Ul+cqb5t9b+IcnbO00/V1Ub5pg7AAAAC6gvt6d5WWf5c621y+fo/44kO6aWVyV56XSdps5NfUmn6S+6e1Jn8LbO8uEZ7GUFAABgiRz0oFpVZye5sNP0oblqWmv3J/lMp+nFM3R9TpKN8xz7qiQ3jjA2AAAAi+CgB9UkTx1a//SIdd1+T1nEsYfHAAAAYBH1IaheMLR+3Yh13X5HVtWpc4y9O8lNY4x9blUN39MVAACARdKHoHpmZ3lvHrn9zFy+Ocs407Xd1lrbN8bY65KcNGIdAAAAE+pDUD2ys7yltbZ3xp4H2jy0Pt09T7tjPzCPOY0ydpKkql5XVVdU1RWbNm2ax1MAAAAwnT4E1cM7y9vnUTfcd7owuZhjJ0laa29vrV3cWrt448aNM3UDAABgRH0Iqms7y3vmUTfcd+00fRZzbAAAABZBH4Lqts7y+nnUDffdNk2fxRwbAACARdCHoLq1s7xhHnXDfbcs8dgAAAAsgj4E1Xs6y4dV1Yzngw4ZvhLvvXOMPZ8r9w73vW8etQAAAEygD0H16qH1M0as6/bbl+TaOcY+Zh4huDv2na21B0asAwAAYEJ9CKpfG1p/2oh13X43tdamu6rv8NhPHWPsq0asAQAAYAH0IahekQMvVvTcEesu6SxfNkOfjw+tzzl2Va1P8owRxgYAAGARHPSgOrUn9COdppdW1awXPqqqZyc5u9P0gRnGvjnJ5ztNr6qqmmNKP5jk0LnGBgAAYHEc9KA65Z2d5aOS/Nwc/d/cWb45yUdHHPucJK+YqWNVrUvyxk7TZ1trDv0FAABYQr0Iqq21D+fAw3TfVFXfO13fqvqNJM/v9m2t7Zpl+Hckub6z/rtV9fRpxl2T5O1JLug0/9JccwcAAGBhrZlP56r61SS/Os2mQzrLl1TVjmn6vLu19tpZhn9dks8kOWZqvA9W1XuS/GUGt4c5K8lrkjynU/PBJO+ebc6ttd1V9ZoM9rquS7IxySeq6p1J/i6De6Sel+T1SZ7cKX1ba234HFcAAAAW2byC6lT/dXP0qRn6rJ2tqLV2bVW9JMlfZRBWVyf54anHdP4+yctba/vmmE9aa5+qqlcm+aMkG5KsT/LTU4/p/GmSfz3XuAAAACy8Xhz6u19r7VNJLkzyniQ7Z+h2a5I3JPmu1tpD8xj7A0mekuSvk+ydodt1SV7VWntla22mPgAAACyiee1Rba29JclbFmUmjzzHnUleUVVHJbk0yalJjkhyV5JrknymtdbGHPu6JN9XVSdkcHubUzLYu3pHkitba1+Y/CcAAABgEvM99HfJtNY2Z3AY8GKMfXeSP1uMsQEAAJhMrw79BQAAAEEVAACAXhFUAQAA6BVBFQAAgF7p7cWUWDpfv+rqoZY5b037sKrk7BOHfo3GuygzAI9BrSUbNxwycv+nnbEmR69fPdmT7my57sZbx6/fcEIeumumu+jN7bBzX5y77ztv/OdPsmf9xlx1z8Vj159ywnG5+8EJ/g2SnNrOyqqjbx67ft3Gc3P8YadNNIcnXPysierXnfkt2X7M+B+HjzjromzasH2iORz3hAvz0CHjj7Hh1DOzatNkd1Vshx2d1NbxB1j16H1fNx76+LSa5W+1tcEHySTZMMf2aVTbmzMeuuHh9TXHnpjVO6cbaDRnn3BGdh412WfYDUcfl52bvvLw+tEv/+XUIetm7N9aS83yMz7WCapMpLXk0GM7b7Zt36wvKgDQVZX82oufOHL/OzdtyZoJjwf76KevyGF5YOz6bzyQ7LvtlrHr7969Ovdcd9fY9UmyfcuW3Hfbg2PXrzlsTbZvHvl29NPatvmhrNu9Zez6XQ9uzsbjxv8ouu/wY/P4kyb7ZdhVD+T45zx17Pq9+/bliKe+dKI57Nu7L2u/8xlj17d9e7Nq3fgBbTBGy6rV4/+/aPtaatWBn//aVVfNXjTX58U5trdancOe/n89vH76zi0Thb61H/yDZOOOsesHbstJv/3+h9e+Nse/gZA6O4f+AgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANAraw72BACAx7bW2sh997WWffsme76v3n107r9j9dj1Z59/bDY/cMbY9SeddXRu237E2PVJcuELz8pRZx01dv3ZTzoht39x3URzeNx5J+SBm+4cu37tukMmev5VW+9LcvxEYwD9JagysSu/9rWDPQUAlrHrr/n6oo5fq1blnHPOfXh90+duy0MPbB97vE179mbzl+8au37n7Vtzx7X3jl2fJEeeeWSe8i3jHxh3xHF78orf+Z6J5pC0rD/ikpF7X3Pttdn3qG8ZHj/RDLZMVD1w64RffNxx19bJJ3HvTZOP8Rj0tauuWrjBLvjuBRnm/oWc02OcQ38BgBWtTboLlgXx6JAKMDNBFQAAgF4RVAEAAOgVQRUAAIBeEVQBAADoFUEVAACAXhFUAQAA6BVBFQAAgF4RVAEAAOgVQRUAAIBeEVQBAADoFUEVAACAXhFUAQAA6BVBFQAAgF4RVAEAAOgVQRUAAIBeEVQBAADoFUEVAACAXhFUAQAA6BVBFQAAgF4RVAEAAOgVQRUAAIBemXdQraojq+rSqvo3VfWeqrq2qvZVVZt6XDbGmJd26ufz+MgYz3VxVf1OVX2lqu6rqoeq6vqq+suqenlVrZvvmAAAACycNfPpXFXXJDknSS3OdBZPVR2W5D8led00mx8/9XhJkiur6lWttS8t4fQAAACYMq+gmuTcRZnFo/3NiP2uGKVTVa1N8sEk39lp3p3k60k2Z/BznTjV/qQkn6iq57TWvjziPAAAAFgg8w2q+21J8sUkn596/JskT12oSbXWXrhQY015aw4MqR9I8rOttduSpKpWJfnBJH+Q5MgkRyT5UFVd0FrbssBzAQAAYBbzDaqvzCCYXttaa/sbq+q1CzqrBVRV5yX5qU7Th5K8rDv/1tq+JO+rqluTfCLJ6iSnJvn5JG9awukCAAA85s3rYkqttT9trV3TDXnLwC/lkUC+O8lPzjT/1to/JHl7p+nnqmrDIs8PAACAjhV9e5qpc1Nf0mn6i/2H+87ibZ3lw5N8z4JPDAAAgBmt6KCa5DlJNnbWPzRXQWvtqiQ3dppevNCTAgAAYGYrPagOX+Dp0yPWdfst2EWiAAAAmFsvg2pV/WFVXVNVW6pqR1XdXlWfrqrfqqqnzWOoCzrLu5PcNGLddZ3lc6tq9TyeEwAAgAn0Mqgm+ZEM7m16eJJ1SU5K8qwkv5jk81X1kao6bYRxzuws3zZ1dd9RfLOzvP/5AQAAWAJ9Dap3Jvlskv+T5PIkDw5t/+4kX6qqp88xzpGd5Qfm8fybh9aPmKljVb2uqq6oqis2bdo0j6cAAABgOn0Jqi3JZUlek+Tk1tpJrbVnttae31r7tgwuiHRpkr/v1ByT5ENVdeos4x7eWd4+j/kM950xqLbW3t5au7i1dvHGjRtn6gYAAMCIehFUW2sfb609r7X2rtbaHdNs39da+3iS5yf5j51NJyT5zVmGXttZ3jOPKQ33XTttLwAAABZcL4LqqNrAz2dwSPB+r6yqE2Yo2dZZXj+Ppxruu23aXgAAACy4ZRVUO7p7UVclecEM/bZ2ljfMY/zhvlvmUQsAAMAElmtQ/WQGt5vZ7/wZ+t3TWZ7PlXuH+943j1oAAAAmsCyDamttdw4Mj8fN0PXqzvIxVTXjRZGGnNFZvrO19sA8pgcAAMAElmVQndI9PHemK/p+bWj9qSOO/bTO8lUjzwgAAICJLcugWlWn5cB7pN45Q9ePD60/d4Sx1yd5RqfpsnlNDgAAgIksy6Ca5IeH1j81XafW2s1JPt9pelVV1Rxj/2CSQzvrH5j/9AAAABjXsguqVXVekl/sNN2W5LOzlLyzs3xOklfMMva6JG/sNH22tebQXwAAgCV00INqVf1AVf12VZ06Qt/nZnAP1aM6zb/WWts7S9k7klzfWf/dqnr6NGOvSfL2JBd0mn9prjkBAACwsNbMp3NV/WqSX51m0yGd5Uuqasc0fd7dWnvtNO1HJvmFJD9fVZ9N8okkX83g1jLbprZfkORFSS4dqn1Xa+1/zTbn1truqnpNko8mWZdkY5JPVNU7k/xdBvdIPS/J65M8uVP6ttba8DmuAAAALLJ5BdWp/uvm6FMz9Fk7Qt0zpx5z2ZfkrUl+eYS+aa19qqpemeSPMrha8PokPz31mM6fJvnXo4wNAADAwjroh/4m+XKS9ye5ZYS+O5L8SZJva639Qmttz6hP0lr7QJKnJPnrJDMdKnxdkle11l45x+HEAAAALJJ57VFtrb0lyVsWcgKttS8m+aEkqaqTk1yY5Lgkx2ZwLupDSTZlcD/TL7bWdk/wXNcl+b6qOiHJJUlOyWDv6h1JrmytfWGCHwUAAIAFMN9DfxdVa+32JLcvwfPcneTPFvt5AAAAmL8+HPoLAAAADxNUAQAA6BVBFQAAgF4RVAEAAOiVXl1MCQBgsR157jFZv33smwjkceccm3V3bBu7/oiN63PT2NUAjw2CKgCw4l17zdUPL//4C48+YNvOnXvzX37zqpHHuuv+7bn7c7fO6/l/9g+elXWHPvKx6/tfdfK86hfDjbfdsQCj3LkAYwA8mkN/AYDHtHXrVi/+cxxq3wDAfAiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABAr6w52BMAADjYjjhq/ch9jznh8Ow8fsMizgYAQRUAeExrreX2D183cv/tX7g977jipw5ou+baq7Nv376FnhrAY5ZDfwGAx7SqmngMIRVgYQmqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0ypr5FlTV+iTfkeR5SZ6W5IlJjk+yNsnmJN9M8tkk72utfXKcSVXV+Ul+JMl3JzktyZFJ7kxyTZL3J3lva23LmGNfnORVU/M/JcmhSW5PcmWS9yb589baznHGBgAAYHIjB9WqOjHJf0nyoiSHz9Dt+KnHxUl+pqo+m+THWmtfH/E51iR5c5I3Jlk9tPmMqccLkvxaVb26tfaxecz/sCT/Kcnrptn8+KnHS5JcWVWvaq19adSxAQAAWDjzOfT3tCT/PI8Oqbcm+cckH0ty7dC2b09yeVU9e8TneGeSX80jIbUluSrJJ5Lc0ul3epK/raoXjDJoVa1N8sEcGFJ3J/lKkk8muavT/qQkn6iqJ484ZwAAABbQuOeofjrJa5Oc1lo7rbX27a2172ytnZfk7CTv6fQ9LMlfVdVxsw1YVW/I4HDf/T6R5PzW2oWttee21k5P8l0ZHKabDPYGv7+qzhhhvm9N8p2d9Q8kOau19uTW2iVJTk7yz5I8OLX9iCQfqqojRhgbAACABTSfoLovyV8meWpr7dmttXe01m4d7tRau7G19ooMDrPd75gMDuedVlUdm+RNnaYvJnlBa+2APbSttY8muSTJ1qmmI5P8+myTrqrzkvxUp+lDSV7WWrutM+6+1tr7knxPkr1Tzacm+fnZxgYAAGDhjRxUW2tfaK3903mcu/nGDA4L3u8HZ+n7M0mO6qy/fqYLGrXWbsiB4fSVVXXmLGP/Uh45F3d3kp9srbUZxv6HJG/vNP1cVW2YZWwAAAAW2KLdnqa1tivJhztNp88S+l7WWf5ca+3yOYZ/R5IdU8urkrx0uk5T56a+pNP0F909qTN4W2f58Az2sgIAALBEFvs+qvcNrR853KGqzk5yYafpQ3MN2lq7P8lnOk0vnqHrc5JsnOfYVyW5cYSxAQAAWASLHVTP7CzvS3LvNH2eOrT+6RHH7vZ7ygx9FmLs4TEAAABYRIsWVKvq0Bx42OzlrbU903S9YGj9uhGfotvvyKo6dY6xdye5aYyxz62q4Xu6AgAAsEgWc4/qv8qBF0h69wz9zuws780jt5+ZyzdnGWe6tttaa/vGGHtdkpNGrAMAAGBCixJUq+pJSd7SabohyR/M0L173uqW1treGfoN2zy0Pt09T7tjPzDiuKOOnSSpqtdV1RVVdcWmTZvm8RQAAABMZ8GD6tQ9Uf8iyfqppr1JXj11FeDpHN5Z3j6PpxruO12YXMyxkySttbe31i5urV28cePGmboBAAAwogUNqlPnpf5Vkid0mn+ltfapWcrWdpanO4d1JsN9107TZzHHBgAAYBEsWFCtqkOS/HmS7+g0v6219ttzlG7rLK+fsdejDffdNk2fxRwbAACARbAgQbWq1iZ5f5IXdpr/IIMLKs1la2d5wzyedrjvliUeGwAAgEUwcVCtqjVJ3pPkxZ3m/5nk9a21NsIQ93SWD6uqGc8HHTJ8Jd7p7tHaHXs+V+4d7nvfPGoBAACYwERBder+on+S5KWd5nclee2IITVJrh5aP2PEum6/fUmunWPsY+YRgrtj39lae2DEOgAAACY0dlCdCqnvTvJDneY/TPLj87hfaZJ8bWj9aSPWdfvd1Fqb7qq+w2M/dYyxrxqxBgAAgAUwVlCdCql/lOTlneY/SvJj8wypSXJFDrxY0XNHrLuks3zZDH0+PrQ+59hVtT7JM0YYGwAAgEUw76BaVasyOLz3FZ3mdyd5zRghNVN7Qj/SaXppVc164aOqenaSsztNH5hh7JuTfL7T9Kqqqjmm9INJDp1rbAAAABbHvILqVEj9n0l+uNP8x0lePU5I7XhnZ/moJD83R/83d5ZvTvLREcc+JwcG7ANU1bokb+w0fba15tBfAACAJTRyUJ3aE/n7SX600/wnSX50wpCa1tqHc+Bhum+qqu+dYR6/keT53b6ttV2zDP+OJNd31n+3qp4+zbhrkrw9yQWd5l+aa+4AAAAsrDXz6PuyJD/RWW9JTkzy13MfTfuwX2itfWWGba9L8pkkxyQ5JMkHq+o9Sf4yg9vDnJXkNUme06n5YAaHHc+otba7ql6TwV7XdUk2JvlEVb0zyd9lcI/U85K8PsmTO6Vva60Nn+MKAADAIptPUB0+b7Ry4J7NUfzWTBtaa9dW1UuS/FUGYXV1BocY//AMJX+f5OWj7M1trX2qql6ZwQWfNiRZn+Snpx7T+dMk/3qucQEAAFh4E91HdaG11j6V5MIk70myc4ZutyZ5Q5Lvaq09NI+xP5DkKUn+OsneGbpdl+RVrbVXttZm6gMAAMAiGnmPamvtXRlc7XdRtdbuTPKKqjoqyaVJTk1yRJK7klyT5DOttTbm2Ncl+b6qOiGD29ucksHe1TuSXNla+8LkPwEAAACTmM+hv0uqtbY5g8OAF2Psu5P82WKMDQAAwGR6degvAAAACKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK8IqgAAAPSKoAoAAECvCKoAAAD0iqAKAABArwiqAAAA9IqgCgAAQK+sOdgTAAA42P7d333vvPp//eqrFmkmACT2qAIAANAzgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9IqgCAADQK4IqAAAAvSKoAgAA0CuCKgAAAL0iqAIAANArgioAAAC9Mu+gWlXrq+qfVNW/q6q/rqobq2prVe2sqrur6vKq+t2qes48xry0qtoYj4+MMf+Lq+p3quorVXVfVT1UVddX1V9W1curat18xwQAAGDhrBm1Y1WdmOS/JHlRksNn6Hb81OPiJD9TVZ9N8mOtta9POM+JVdVhSf5TktdNs/nxU4+XJLmyql7VWvvSEk4PAACAKSMH1SSnJfnn07TfmuS2JA8lOSXJuZ1t357k8qp6YWvtU/N4rr8Zsd8Vo3SqqrVJPpjkOzvNu5N8PcnmDOZ84lT7k5J8oqqe01r78ojzAAAAYIHMJ6h2fTrJu5J8pLV2a3dDVZ2V5DeSvHyq6bAkf1VV57XW7h1l8NbaC8ec10zemgND6geS/Gxr7bYkqapVSX4wyR8kOTLJEUk+VFUXtNa2LPBcAAAAmMV8zlHdl+Qvkzy1tfbs1to7hkNqkrTWbmytvSKDw2z3OybJGyea6Ziq6rwkP9Vp+lCSl+0PqUnSWtvXWntfku9Jsneq+dQkP79kEwUAACDJPIJqa+0LrbV/Oo9zN9+YwWHB+/3gfCa2gH4pj+w53p3kJ1trbbqOrbV/SPL2TtPPVdWGRZ4fAAAAHYt2e5rW2q4kH+40nb7UoW/q3NSXdJr+orsndQZv6ywfnsFeVgAAAJbIYt9H9b6h9SMX+fmGPSfJxs76h+YqaK1dleTGTtOLF3pSAAAAzGyxg+qZneV9SUa6mNICeurQ+qdHrOv2Gx4DAACARbRoQbWqDs2Bh81e3lrbM2LtH1bVNVW1pap2VNXtVfXpqvqtqnraPKZxQWd5d5KbRqy7rrN8blWtnsdzAgAAMIHF3KP6r5Ic1Vl/9zxqfySDe5senmRdkpOSPCvJLyb5fFV9pKpOG2GcMzvLt7XW9o34/N/sLO9/fgAAAJbAuPdRnVVVPSnJWzpNN2Rwj9JR3ZnB3s9tGZzXel4OPL/1u5N8qape2Fq7fJZxujUPzOP5Nw+tHzFTx6p6XZLXTa1uveiii66Zx/MAAACsJGcsxCALHlSr6tgkf5Fk/VTT3iSvnroK8ExaksuS/GGSv2mt3TE05qoMLoz0piTfOdV8TJIPVdW3Tnc/1ymHd5a3z+PHGO47Y1Btrb09B97SBgAAgAks6KG/U+el/lWSJ3Saf6W19qnZ6lprH2+tPa+19q7hkDq1fV9r7eNJnp/kP3Y2nZDkN2cZem1neaTzY2fou3baXgAAACy4BQuqVXVIkj9P8h2d5re11n57oZ6jDfx8kv/TaX5lVZ0wQ8m2zvL6GfpMZ7jvtml7AQAAsOAWJKhW1dok70/ywk7zH2RwQaXF0N2LuirJC2bot7WzvGEe4w/33TKPWgAAACYwcVCtqjVJ3pPkxZ3m/5nk9a21Nun4M/hkBreb2e/8Gfrd01mez5V7h/veN49aAAAAJjBRUJ26v+ifJHlpp/ldSV67iCE1rbXdOTA8HjdD16s7y8dU1YwXRRrSvVLVna21B+YxPQAAACYwdlCdCqnvTvJDneY/TPLj87hf6SS6h+fOdEXfrw2tP3XEsZ/WWb5q5BkBAAAwsbGC6lRI/aMkL+80/1GSH1uKkFpVp+XAe6TeOUPXjw+tP3eEsdcneUan6bJ5TQ4AAICJzDuoTt3T9F1JXtFpfneS1yzRntQk+eGh9Wlvf9NauznJ5ztNr6qqmmPsH0xyaGf9A/OfHgAAAOOaV1CdCqn/MwcGxT9O8uqlCqlVdV6SX+w03Zbks7OUvLOzfE4ODNjDY69L8sZO02dbaw79BQAAWEIjB9WpPZG/n+RHO81/kuRHJwmpVfUDVfXbVXXqCH2fm8E9VI/qNP9aa23vLGXvSHJ9Z/13q+rp04y9Jsnbk1zQaf6lueYEAADAwqpRL85bVT+U5L2dppZBaJwtJA77hdbaV4bGfXWS/zU13meTfCLJVzO4tcy2DM5FvSDJi5JcOjTeu1prrxlh7s9O8tEk66aadmSwp/XvMrhH6nlJXp/kyZ2yt7XW/uXoPxoAAAALYT5B9dUZBMpJPK+1dtkCjLsvyVuT/HJrbc8oBVX10gwu+LRhrr5J/jTJj8yxpxYAAIBFMNF9VBfIl5O8P8ktI/TdkcHhxt/WWvuFUUNqkrTWPpDkKUn+OjPvBb4uyataa68UUgEAAA6OkfeoLoWqOjnJhUmOS3JsBueiPpRkUwb3M/1ia233AjzPCUkuSXJKkvVJ7khyZWvtC5OODQAAwGR6FVQBAACgD4f+AgAAwMMEVQAAAHpFUF0BqurPqqpNPe6YR93LOnWtqu6eul/uKLX/dKj2+8b/CQAAYOWoqrcMfVYe53Hmwf45DiZBdWX4WGf5cVV1/oh1lw6tH5/BxazmW7s3ySdHrDuoqurS7gvAwZ4PAACsFFV15kKF7TULOC8OnsuG1i9NcvUIdZfO0HblPGs/31p7cIQaAAB4LPqbMWq2L/gslhFBdQVorX2tqu5OcsJU06VJfm+2mqo6PskTp1bvy+B2QPtr3zZH7TFJLuo0fWymvgAA8FjXWnvhwZ7DcuPQ35Xjss7yc0fof2mS/eej/nGSu/bXjnCe6nM7tYmgCgAALCBBdeW4rLP8uKp64kwdp1w6VLv/HNPjkjxpHrW7k3xqztkBAACMSFBdOYb3al46R//921sGIfXjY9QmyeWttW3TdaqqQ6vqx6rq/VV1fVVtrqrtVXVzVf3vqvoXVbVhjufqjndUVf1UVf1/VXVTVW2tqj1VtaWqbqyqv6+q/1hV31tVa4dq3zV18aSPDbXPdJW1d406LwAAWCpVdXRV/bOq+u9V9emqurOqdkx9zr5t6jPxm6vq1HmOu76qXjV1R5Hrpz5j762qbVOf3z9ZVb9bVT84/Bl+/1WOk9w4NOyNM3zWvmyu+ThHdYVorV1dVXcmedxU06VJ/sd0fYfOT/1qa+2+qhoOqr87Q+1I56dW1SuT/PskJ0+z+bSpx/cm+ZWqel1r7X9PN05nvO9P8s4Mrkw87PCpx5lJnpfk3yT5D0l+YbYxAQBgOamqf5vkF5McMkOXk6cez8vgc/b/k+QtrbVZ73ZRVc9M8idJzppm84apx2lJnp3kZ5K8P8kPjfMzjEpQXVkuS/LPp5ZnO0/10jxyjun+gHplHrmo0iVVVTP8Qs95fmpV/UaSXx5qvivJDUl2ZRAoz5xqPznJX1XVj7XW/mi6yVbVpUn+PAf+vt6f5LokW5IcmuTEJGfnkaMEho8W+GoGV1s7JsnTO+0zXYHtqzO0AwDAwXJBDgypdya5OYPPxOsz+Dx80tS2tUnelMGOrNfPNODUrS3/NoMdP/ttSXJNkgeSrMvg9MAnTI2ZPPqz9vUZfK4+NMklnfZPZPqrF39lpvnsJ6iuLB/LI0H1xKq6oLV21TT9Lu0sX5YkrbVWVZ9K8pIMfhEvyvS/QN3aXUn+obuxqn4yB4bUv03yq621y4f6PSWDvbbPTrI6ye9X1Rdba9MFxP+cR35Xr0vyk0k+Nhykpw5B+M4kr8jg3NmHtdbemuStU6H3Y512V2ADAGC52JvkQ0nem+RvWmv3DHeoqouSvCXJD0w1va6q/ndr7YMzjPn/5JGQeleSf5Hkg621PUPjHpLkOUlelkEgfVhr7Y+T/PHUfVO7h//+aGvtplF/uC5BdWW5bGj90iSzBdWWwbcc+308g6C6v89cQfWzrbWHvyGpqjMyCJX7/bck/3K6PbOttS9V1Xcm+UgG4XJ9kt9K8n3dflPH1j+lM98XtdaunWZeaa09lMEf7oeqavV0fQAAYBn78dba1tk6TO34eWlV/V4e2ZP6i0keFVSrak2S7+k0vaq19nczjLsryf9J8n+W4rO2iymtIFMB7rZO06XDfabOT71gavWq1tq9nc2zXlBphPNTfy6DwJkMDp392dmOh2+t7U7y40n2f1vzPVX1+KFup3WW754ppE4z9t5R+gEAwGKb5QKeMz3+y3TjzBVSh/xCkh1Ty8+qqhOn6XN8Bof27vfJafpMN49F/6wtqK483bA53Xmql3aWLxva9qUkm6eWL5nmfqoznp9aVauS/Ehn238e5Rd46lCA/Xt1K8k/GerSPab9hKo6Za4xAQDgsa619mCSr3Wavm2absPnjz518WY0Pw79XXk+lsE5mskg2A2fp3ppZ7kbatNa21dVn87garzH5tHnqXZrdyT5bGf9oiQbO+sfncecv5LB4b9JcnGSt3e2XZVkW5LDMgiyH5q6SvDlAQCA5WGmC3jO5Otzdaiq8zLYkfSkDPaMHpFH57uzO8uPul1Na+2Bqro+gwslJcmfVtVrW2vz+Sy/KATVlWf4cNzn5cDzVC/tLH88j/bxDILq/tqZgupnWms7O+vfMjTOHzx6h+yMntBZPuD2M621XVX1O3nkAk1PSfK5qro2gz/4T0/N5eZRnwwAAJbSQl7As6qeneQ/JnnGPEuPnqH93+eRHUVnJvm7qrolg2vJfDKD69JcN/+ZTkZQXWFaazdM/WLtP7fz0gwuajR8furXW2t3TzPE8HmqvzNVO9f5qccOrX/3fOc+5ahp2t6c5JQkP9ppO3fq8S+n5veNDG5h847W2jVjPjcAAPRWVb0uye/lwNPxRrVuusbW2h9U1VlJfqkz7mlJXjv1SFXdnsHFmN7ZWrtijOeeN+eorkyXdZaf2znX9NIZ+nR9Psn+k7S756nOdf/Uw+Y9y+k96neytbantfbqDA4P/qskO4f7ZHBYw/+d5Kqq+m9VNe0fIgAALEdV9dQk/yOPfCa/P8lbMzga8pwkRyY5pLVW+x+Z/gjKR2mt/XIGp+D9SR7JAl0nZ3CLyMur6n1VdfQkP8so7FFdmT6W5FVTy/v3on4tcx/2m9banqr6TJLvSnJMBof0fnmo9qEknxsqfaCzvLm1dvRYM59Fa+1jST5WVYcmeWaS78jghsLPziNXG16Vwb2fjkvyzxZ6DgAAcJD8Sh7ZqXNTku9ord0+R80Row7eWvtCkh+uqrVJnp7BZ+3nZJADuuO8LMmpVfWcxbz6rz2qK9NlQ+vPm/rvpZ222b5d6W6brvbTU/dR6rqzs3xUVZ0w+xTH11rb3lr7+9bar7fWviuDUPqaJLd2uv1QVT1rseYAAABLZeoox+6pdf92hJCaTHMBpbm01na31v6htfYfWmsvzuCz9ktz4AWenplF3ikkqK5ArbUbk3yz03Tp0Pmp17TW7nx05cMOOE91hPNTk+QfhtZfMOp8J9Va29Zae9fUc+7ubJruPNl93ZVpbsEDAAB9c0ySwzvrw0c3PkpVnZtk4p1HrbVdrbU/z2AH1n2dTXN+1s5459ImEVRXsss6y8/NI3tGh7dN53N55J5Kl0zVznZ+alprd2RwiPB+PzPiPBdMa+3rOfCbnsdN023b0PqhizcjAABYEGvHqPnxhZxAa+2uDO64sd+iftYWVFeubpg8LslPd9ZnPal66rDef5xa3ZjkX3U2b00y05W+/kNn+RlV9a9Hmuksxtjj2f2m6f5ptt8xtH7OPMcHAICldm+S7ql3l8zWuaouyIGf4Wfqt9CftR9IsqOzPvZnbUF15Rre69n9ZR7l6l/dPt3aT7bW9sxQ854c+C3LW6vqV6pq1ot2VdWhVfXyqvr8NJtfWVXvrapnzjXhqvoXOfCmxtPt+b09B55P+6+ravVcYwMAwMEy9fn7k52mN03dUuZRquqiJH+TRy42OptLqurDVfWCuT4TV9X358Dr1kz3WXtvki91mv5FVY0yj0dx1d8VqrV2c1XdmGT4F/i6EU+8ninMTnd+6v7n3FdVL81gb+wZGXwR8u+SvL6q3pPk8iT3TLVvTHJeBlcUe0Fmvr3NmiQ/lMHFkb6ZwY2Hv5DkliQPZnA4wblJ/mmS53fq/iHJR2cY848zuJVNkrw6yfdW1Vdz4KW4/7619l9n+lkBAGBUVfWRMcp+r7X2l531tyb5J1PLj0vyhar6vSSfyOBz7MkZ3KrmFRl8hv5iBtdv+bbZppbkhVOPu6bmeUWSG5NszuCQ47OTfF+Sl+SRHZ03ZPCZejp/nOTbp5ZfkOSOqvpSBp/d21T7la21X51lXoLqCvexPDqoXjZi7WczOLzgkGnGnFFr7a6qekaS9+WRPbGnJfmFEZ93Nmckef0I/b6c5Adba8Mnc+/3bzP4I3/q1PoJeeSPfr8HxpkgAABMY7oLD83lgHDbWvtwVf2nJG+Yajo6yS9NPYZ9I8kPJHnXPJ7vxCQ/OvWYzc1JXtRae2iG7b+f5EUZhN/987x0qM/Rc03Gob8r23ShctSb/m7Po68m9mAG38zMVXtXBhdg+ucZfCPTZq/INRl8Q/Tkabb9fZLfyiB8zhQ897s5yRuTPGPq4k4zzW9LBpfUfn2SD2dwW5vtM/UHAIA+aK39myT/IsndM3TZluTtSZ7aWrtphCG/muRNGeyk2j1H33uS/PskF7XWrp5ljnsy2AP7iiR/kcHe2W2ZOxMcoFqbV3+Yt6lb43xHkpMyOOR3TwZ7LL+RwW7/2W6V0x3niCRPyeDwg+MzOOz3oSR3ZRBkr2x+oQEAWOGqal2SZye5MIMLHN2Xwalxl82yp3OuMQ/NYMfREzI44nBDBhdGujfJlUm+NMu1ahacoAoAAECvOPQXAACAXhFUAQAA6BVBFQAAgF4RVAEAAOgVQRUAAIBeEVQBAADoFUEVAACAXhFUAQAA6BVBFQAAgF4RVAEAAOgVQRUAAIBe+f8DcONBCWETLI0AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots(1, 1, figsize=(15, 10))\n", "jmin = 210\n", "gridding = 2\n", "sisal = 30\n", "\n", "v2max = 0.01 / otransports[osals <= sisal].sum() / m2lon * 24 * 365 * 4\n", "print ('oo', otransports.sum()/ (24 * 365 * 4), otransports[osals < sisal].sum()/ (24 * 365 * 4), \n", " otransports[osals >= sisal].sum()/ (24 * 365 * 4))\n", "c, xedge, yedge, im2 = ax.hist2d(olons[osals < sisal], \n", " odepths[osals < sisal],\n", " weights=otransports[osals < sisal],\n", " bins=[np.array(f_lons[jmin:314:gridding]), new_ws], cmap='Greens_alpha', \n", " vmax=v2max, density=True,)\n", "\n", "print (c.max()/v2max)\n", "\n", "vmax = 0.01 / transports[sals < 31.2].sum() / m2lon * 24 * 365 * 4\n", "print ('oo', transports.sum()/ (24 * 365 * 4), transports[sals < 31.2].sum()/ (24 * 365 * 4), \n", " transports[sals >= 31.2].sum()/ (24 * 365 * 4))\n", "\n", "c, xedge, yedge, im = ax.hist2d(lons[sals < 31.2], \n", " depths[sals < 31.2],\n", " weights=transports[sals < 31.2],\n", " bins=[np.array(f_lons[jmin:314:gridding]), new_ws], cmap='Oranges_alpha', \n", " vmax=vmax, density=True)\n", "print (c.max()/vmax)\n", "\n", "vmax = 0.01 / transports[sals >= 31.2].sum() / m2lon * 24 * 365 * 4\n", "c, xedge, yedge, im = ax.hist2d(lons[sals >= 31.2], \n", " depths[sals >= 31.2],\n", " weights=transports[sals >= 31.2],\n", " bins=[np.array(f_lons[jmin:314:gridding]), new_ws], cmap='Purples_alpha', \n", " vmax=vmax, density=True,)\n", "\n", "print (c.max()/vmax)\n", "\n", "v2max = 0.01 / otransports[osals > sisal].sum() / m2lon * 24 * 365 * 4\n", "c, xedge, yedge, im2 = ax.hist2d(olons[osals >= sisal], \n", " odepths[osals >= sisal],\n", " weights=otransports[osals >= sisal],\n", " bins=[np.array(f_lons[jmin:314:gridding]), new_ws], cmap='Blues_alpha', \n", " vmax=v2max, density=True)\n", "\n", "print (c.max()/v2max)\n", "\n", "\n", "\n", "ax.invert_yaxis()\n", "ax.set_ylim(250, 0)\n", "ax.set_xticks([])\n", "ax.text(-123.57, 270, \"West\")\n", "ax.text(-123.1, 270, \"East\")\n", "pc = draw_patches(fmask, jmin, 314, np.array(u_lons), np.array(w_depths[0]))\n", "ax.add_collection(pc);\n", "print ((mymesh.e1f[0, igrid, jmin:314].sum()).values)\n", "fig.savefig('oroberts_final.png')" ] }, { "cell_type": "code", "execution_count": 21, "id": "b24a0c90-fa35-4e94-9405-73e3e85a9d69", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "15\n", "16\n", "17\n", "18\n" ] } ], "source": [ "month = {1:'jan', 2: 'feb', 3: 'mar', 4: 'apr', 5: 'may', 6: 'jun', 7: 'jul', \n", " 8: 'aug', 9: 'sep', 10: 'oct', 11: 'nov', 12: 'dec'}\n", "lons, depths, transports, sals = np.array([]), np.array([]), np.array([]), np.array([])\n", "months = 4*np.ones_like(lons)\n", "years = 15*np.ones_like(lons)\n", "for y in [15, 16, 17, 18]:\n", " print(y)\n", " for m in month:\n", " newlons, newdepths, newtransports, newsals = get_data_backward_vic(month[m], y)\n", " newmonths = m*np.ones_like(newlons)\n", " newyears = y*np.ones_like(newlons)\n", " lons = np.concatenate((lons, newlons))\n", " depths = np.concatenate((depths, newdepths))\n", " transports = np.concatenate((transports, newtransports))\n", " sals = np.concatenate((sals, newsals))\n", " months = np.concatenate((months, newmonths))\n", " years = np.concatenate((years, newyears))\n", " \n", " newlons, newdepths, newtransports, newsals = get_data_backward_gi_vic(month[m], y)\n", " newmonths = m*np.ones_like(newlons)\n", " newyears = y*np.ones_like(newlons)\n", " lons = np.concatenate((lons, newlons))\n", " depths = np.concatenate((depths, newdepths))\n", " transports = np.concatenate((transports, newtransports))\n", " sals = np.concatenate((sals, newsals))\n", " months = np.concatenate((months, newmonths))\n", " years = np.concatenate((years, newyears))" ] }, { "cell_type": "code", "execution_count": 22, "id": "8d0b7c36-b344-42fd-9ef6-26602b3a7efb", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "15\n", "16\n", "17\n", "18\n" ] } ], "source": [ "month = {1:'jan', 2: 'feb', 3: 'mar', 4: 'apr', 5: 'may', 6: 'jun', 7: 'jul', \n", " 8: 'aug', 9: 'sep', 10: 'oct', 11: 'nov', 12: 'dec'}\n", "olons, odepths, otransports, osals = np.array([]), np.array([]), np.array([]), np.array([])\n", "omonths = 4*np.ones_like(olons)\n", "oyears = 15*np.ones_like(olons)\n", "for y in [15, 16, 17, 18]:\n", " print(y)\n", " for m in month:\n", " newlons, newdepths, newtransports, newsals = get_data_forward_vic(month[m], y)\n", " newmonths = m*np.ones_like(newlons)\n", " newyears = y*np.ones_like(newlons)\n", " olons = np.concatenate((olons, newlons))\n", " odepths = np.concatenate((odepths, newdepths))\n", " otransports = np.concatenate((otransports, newtransports))\n", " osals = np.concatenate((osals, newsals))\n", " omonths = np.concatenate((omonths, newmonths))\n", " oyears = np.concatenate((oyears, newyears))\n", " \n", " newlons, newdepths, newtransports, newsals = get_data_forward_gi_vic(month[m], y)\n", " newmonths = m*np.ones_like(newlons)\n", " newyears = y*np.ones_like(newlons)\n", " olons = np.concatenate((olons, newlons))\n", " odepths = np.concatenate((odepths, newdepths))\n", " otransports = np.concatenate((otransports, newtransports))\n", " osals = np.concatenate((osals, newsals))\n", " omonths = np.concatenate((omonths, newmonths))\n", " oyears = np.concatenate((oyears, newyears))" ] }, { "cell_type": "code", "execution_count": 24, "id": "cb133f81-94a6-46f5-b2f4-a4a452d5d6a3", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "oo 1112208465.1924486 935953646.5997791 176254818.5926676\n", "5.446429832539066\n", "oo 906491677.6605598 681792440.4049758 224699237.2555846\n", "2.45928290351305\n", "0.9165026140189085\n", "0.8639584274433374\n", "35823.448754892044\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAwAAAAGUCAYAAACC6UXyAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAAzFUlEQVR4nO3deXhcZ3n38d89Gq22JO+WvO927DiGEAhblpaQAqFAISmFkLxhS+AFyr60DUvL0tIWWgpvrzaQQoCSBkihaQqENQmhgSwkcRxvcRwvsSPvi2RZ1jL3+8eM8fFYI52RnpGOdL6f65rL5zzznHseydLo+c3ZzN0FAAAAIB0yoz0AAAAAACOHAAAAAACkCAEAAAAASBECAAAAAJAiBAAAAAAgRQgAAAAAQIoQAAAAAIAUIQAAAAAAKZK6AGBm88zsL8zs12a228xOmNlOM7vLzN5tZtNHe4wAAABApVia7gRsZu+S9FlJ9QN0Oyjp7e7+7ZEZFQAAADByUhMAzOyTkq4van5c0m5JcyQtLnruTe7+1ZEYGwAAADBSUhEAzOzVkm6NNK2XdJW7/zbS5zxJX5d0VqGpV9IL3P2+ERsoAAAAUGHjPgCYWbWkjZIWFZqeknSOux/qp+8USWslzS40/dLdLxyRgQIAAAAjIA0nAb9epyb/kvS+/ib/kuTuByW9L9J0gZkRAAAAADBupCEAXBFZ3i3pe4P0/89Cv/62BwAAAMa0cR0AzKxe0iWRph+5e+9A2xSevyPS9IpKjA0AAAAYDeM6AEhaKak2sv6rmNtF+80rnBsAAAAAjHlpCABRj8fcrrhfcR0AAABgTBrvAWBB0fqOmNttH6QOAAAAMCZlR3sAFdZUtH445nZHitYbS3U0s2slXStJ9fX1z1q4cGHsweXctX3frjPaJ9Q3qru3J3adgTTUTVBX9/EwtWrD1aqvnaAT3Z1hatXUqztUrdo6He86FqRWTe0EdXSFGZdVTVTnia4gtaRanegKVSsjeV+QSqaMZIEuS2wZZaoClLGsMpkwv4tVVTWSh/m+Z7N16u0N87uYra5Xb0+YWtXV9erpCfMzX13doJ6g4wpXK9z3K9zXmK2uD/YzUVNdF2xcdTX16ukJ857aUNsQ7L1+Ql29jgd6f26sq9eJ3jC/2w3V1cpkA7x5Sap2V1VVoGleLqdMdV2QUlZdH6SO3GUh3ugLtWQWplZEV1eXhnu5/fXr1+939+mBhvQ74z0ATCxaj/uOVtyvZABw9xsk3SBJq1at8ltuuSX+6CS9/ksfOqPteat/Tzv7CQZDsWbp+Xp852NBap2z5DnavH1tkFpnLzlPW3c8GqbW/LO1fccDQWqds3ClHtv8myC15i44T/dsDvP9qp30XD3w+OYgtdwWaPvmLUFqyRvkJ/q9qm7ZMpl6KRsmTGSzE1XbNOD5/rHU1k1XQ8PeACOSGpvnyrs3Bqk1bfoKHdgX5vdnxsxV2rsnTK2WlrO1py3Mz3xr6znavTvMuFpbV2tPW5has2ado7anw3yNs2Y9U08/vT5IrZktq7Rv77ogtWa1rtT+PY8EqbVk7irt3f3bwTvGsGbxGm3f+VCQWucvPUcbtoQZ10VLVmrLgTC/2+e2zNaE6cWfXQ7NjO4eTWoOM2+sOtGp2tazBu84mL5e1bacPfw6ktxzqqopnuYNuZgyVYHChPS7MLFu3fB/J1evXl18VEoQ4/0QoOqi9bgzguKP/IrrDChu2hvvN2EDAABA8oz3AFC83zHu/qvi/VNl7b+0mLuR4vYDAAAAQhnvAaCjaL0h5nbF/doDjAUAAAAYdeM9AOwrWm+NuV1xv/3lvCiHAAEAACCpxnsAKD4rZ37M7Yr7lXV2D4cAAQAAIKnGewAovvzNuTG3i/brllTWJVPYAwAAAICkGtcBwN13Stoaaboo5qbRfve4l3ehc/YAAAAAIKnG+30AJOl7kt5fWL7YzOa5e8k7ApvZPJ0eAG4t9wXdPdbknj0AAAAAY8/6DRuUy+VGexhDNq73ABR8VdLJ/6GMpI8O0v9jOvV96ZD07XJfkD0AAAAA49dYnvxLKQgA7v6YpG9Gmt5iZm/pr6+ZXSfpzZGmv3f3sq4AVHjNoP0AAACAUNJwCJAkfVDSBZIWFta/bGZ/KOk/JO2WNFvS6yS9PLLN/ZL+bigvxh4AAAAAJFUqAoC77zWzl0m6Q9K8QvMrCo/+rJX0cnfvHOLrcQ4AAAAAEmncHwJ0krtvlLRa0j/rzDsEn3RA0qckPdvd9w71tdgDAAAAgKRKxR6Ak9z9qKR3mNkHJF2s/A2/Jit/p98nJN3t7r0BXoc9AAAAAEikVAWAk9z9uKQfVqo+ewAAAACQVKk5BAgAAABASvcAVBqHAAEAAIw9Y/0GX3GxB6ACOAQIAABg7EnD5F8iAFQENwIDAABAUhEAKoA9AAAAAEgqAkAFsAcAAAAASUUAqAD2AAAAACCpCAAVwB4AAAAAJBUBoALYAwAAAICk4j4AFcB9AAAAAJJl46ZNqbnM52DYA1AB7AEAAABIFib/pxAAKoBzAAAAAJBUBIAKYA8AAAAAkooAAAAAAKQIAaACOAQIAAAAScVVgCqgnEOAXnDeM89oX7P4Ward2hRkLEvnLFXGe4LUWjF7gWrUFaTWqlkLNMFOBKm1fO4CTbL2ILUWz1qimp4jQWpNm7VER7rDfO994jx19IQ5ZCxbO1011WGyf+PEJh05dChIrRnTp+hYR0eQWrNaJqvr+LFh15k9Y7KOdU4a/oAkzZ7WogMHw5yAtmj2Im1v7A1Sa8mcpdraGGZcy+cu09bGMB9srJh/lrZNrgpSa8ns5doyOUgprVl4tva1hPnTuWDWCu2eHqSUprUu1RO7wvw/Lp6zTAdaw7x3LZmxWD2Tw7xHTJ/SqiOHnwxSq7luqubaxCC1pniTZreH+X1sbqxVw9HHg9Sqb5yh2prhvw9KkmWyMg3/b5CrSlKgDz9zOYkPUoeMAFABLo/1i+Jyvf1Vrz2j/as/uEOPbN4YZCy57k5t2r42SC3r7dCm7Y8FqVWb69S6rY8EqdXQd1gPb34wSK2q7gO6b1OYWmcd79TaTeuD1Fq4OKv2Q5uD1Jo689lqqHoqSK2a7EJVV4epla1qUE12R5BatZlaeW74v0O1WqoDh+4LMCKpt2G5du+8N0itSdl27dr5cJBaU2uPqy1QrZl1PWrbHabWvEbXwbYw7xEnJvRq3677g9Sqn2qyfWHeU2c29al9X5ifr9kNHdoY6Ps1o/Gwek5sD1Jr1sE9Mt8fpFbd9j3q27ctSK1cV4f6ntoapFbH8WOanQnzwZE/+ZCmLp4dplbHPjWsel6QWsr1KjNj5bDLuOeUqakPMCBJLlkm/2HW+g0buMJPmTgEqALipuQQaRoAACDNmPyXjwBQAR5z91bcfgAAAEAoBIAKYA8AAAAAkooAAAAAAKQIAaACOAQIAAAASUUAqAAOAQIAAEBSEQAAAACAFOE+ABVQzn0AAAAA0L8NGzdxmc8KYA9ABXAIEAAAwPAx+a8MAgAAAACQIgSACuAqQAAAAEgqAkAFcAgQAAAAkooAAAAAAKQIAQAAAABIEQIAAAAAkCLcB6ACuA8AAABAaVzff3SxB6ACOAkYAACgNCb/o4sAAAAAAKQIAQAAAABIEQIAAAAAkCIEAAAAACBFCAAAAABAihAAAAAAgBThPgAVUM59ADJtZ2awN5/7Ur353JdKkrp6uvX+//h/wccIAABQCRu3PKmcc6+jJGMPQAWEvA9AXXXNcIcDAAAwYpj8Jx8BAAAAAEgRAgAAAACQIgQAAAAAIEUIAAAAAECKEAAAAACAFCEAAAAAAClCAAAAAABShAAAAAAApAh3Ah4D5i5eMORtly9eoZopYcZxzqLVaphaFaTWqvkrVd0U5kYhq+asVM2EviC1zm5ZqgmNuSC15s04R+11Yb5fy2cv17T62iC1ZsxcqK3N1UFqTZs6R9v31geptbR1jo4cClNrYctsHd8//P/HmTMX6NCkrgAjklpal+qshjA/pysXrNb2SWF+tpbNOkcHq3qC1GpoWa7VR8J8vxrrF2rpxKNBak2um6XzlswPUssyrertOhik1v5DLartXhKk1sGOeVp+Yl+QWl3H5qpmy64gtXYvWaLDG/YHqTVh6Xk66GHeu5qmrFFXe5j31Pr5K1TT9osgtapbFsgmhLkBqGUyksW7MengxQLVQWIQABLO3bXz8NYhb9+9uUNz50wPMpaH9vxWzY11QWqt3fuwujNhJh2P7n1EXVXdQWo9cWiTqpobgtTacugxbT809P+7qNa6Bh19elOQWsuaatR89LEgtSY3SzvbHw1Sq2l6lboPPxykVnNTj/ra7h92ncaa49q168EAI5Ky1Se0d/faILWWTczowK4w3/fq6ox6dj4cpFbW69W1Y32QWjVV9ercHubntKb3hOY07gxSq+NQlY5s3RCkVm+uRh1bw/xM9Paajm0LMy5ZjdqfDFPLM3U6vH1jkFo91ZN1ZMfmILVyVqPjW9cFqdVzvF0rZ4SZINuRnWq49P8GqeVyZWauDFJLkixTpfXr1yuXC/MhGUYXhwAlnJG6AQBAAjD5Hz8IAAAAAECKEAAAAACAFCEAAAAAAClCAAAAAABShAAAAAAApAgBAAAAAEgRAgAAAACQIgQAAAAAIEUIAAAAAECKEAAAAACAFCEAAAAAAClCAAAAAABShAAAAAAApAgBAAAAAEgRAgAAAACQItnRHgAG5u6jPQQAADBGbdi4SblcbrSHgYRhD0DCmdloDwEAAIxRTP7RHwIAAAAAkCIEAAAAACBFCAAAAABAihAAAAAAgBQhAAAAAAApQgAAAAAAUoQAAAAAAKQIAQAAAABIEQIAAAAAkCIEAAAAACBFEhMAzGyBmfkQHhvLfJ0VZvYZM3vQzPaaWZeZbTOzO8zsLWbWWKmvEQAAABht2dEewEgxs6ykj0v6M0lVRU/PLzwulfRRM7vG3X8xwkMEAAAAKi7JAeBuScdj9NsZs96Nkq6OrLukDZL2S1ooaW6hfZ6kH5vZZe7+45i1AQAAgDEhyQHg/7j7thCFzOx9On3yf7ekt7r75kifSyTdJGmW8t+X75jZOe6+PcQYAAAAgCRIzDkAlWJmUyV9LNL0kKRLo5N/SXL3n0q6UFJHoalJ0idHZJAAAADACBn3AUDSOyU1R9avc/cT/XV09yd0+qT/SjNbUMGxAQAAACMqDQHgisjyfe5+/yD9vyKpq7CckfSaiowKAAAAGAXjOgCY2SJJqyJNtw+2jbsflHRvpOkVoccFAAAAjJZxHQAkPbNo/Vcxt4v2e0aYoQAAAACjL8kB4G/N7FEzO2xmJ8yszczuN7N/NLOLYtZYWbT+eMztov2azGxOzO0AAACAREtyALhC0tnKn8BbI2mmpPMkvVvSnWb2GzM7e5AaCyLLfZJ2x3zt4kt/LuivEwAAADDWJDkAHJB0n6SfSfqN8jfsinqOpPvM7OUD1GiKLLe7e1/M1z5StN4YczsAAAAg0ZIWAB6U9A5Ji9x9mruf7+6XuPtz3X268nsAvhvpXy/pFjN7Rol6EyPLce4qXKpvyQBgZtea2QNm9sChQ4fKeAkAAABg5CXmTsCFu/6eN0ifByVdYWbvkvRPheaGwvKF/WxSHVnuLWM4xX2r++2VH9MNkm6QpFWrVnkZrxHb0pb5Q972+UtWqzvTNXjHGOY1zdLhE2FCztwJ07S3Y1+QWksnztSeY21halU3qbNzT5Bazc3T1DK7O0it+VNatby7Y/COMdRNblXjkeKdXEMzYdosTe85GqRWw4SZqq1bFKTWgX2TdfDwimHXaWtvVnXnsgAjkh7pnKZZnQuD1Nqza6bm9cwOUuupzgU6mOn31ihl2/50i450rQ5SS9ubNbNpTZBSR+vnafNTYd67ulvO0rGmuiC1Dra3qC/7jCC1dj89TVkb8E9obAcOtqhlzvlBamn6WZrYG+ZPY8O8NerpqwpSq2nuMs2pCvM3qHbGXGWq4h5hPDCryqoiEwmgSGICQDnc/YuFT/3fVGi6wMzOc/cHiroeiyyX845d3PdYv71GyHvOvXzQPmamhUsXnNHedrxNGQuzo6e966hW1Q09jER1HT+qhTPDnFvd1bFPy5qWBKmVPdymlqlhJlY62KZzWucGKXVs13rVde0NUqtm613q3PhokFrNfZ06sTtMraYTferc9HCQWrkZ56tz15Zh1/FJpu5924Y/IElV02rUuz9MrZxVq6ttbZBaNblGHdm6Lkit3mkN6tzzZJBa2Smuo3sfCVKr73iHcjvC1Mp4izp3xb2exMByU6vVvXdrmFpNWfUeLD59bWiq+6Sj3WF+JnI5V/eOh8LUUrWOPflwkFrWtV8t2hGm1s7dav3Yl4PUcndlsrWx+2/YuFG5XK7UyKQDG4KMC+NP0g4BKseni9Zf0k+f6EemDWXULu7bXsa2o8KdzwwAAEiT0pN/YGBjNgC4+1adfrWe/vb3R/fvTTCzuCfzthatF5+ADAAAAIxJYzYAFDwdWZ7Wz/Mbi9bjHr8S7ZeTtLmcQQEAAABJNdYDQPRQnf6u8vNY0fq5MetG+21z93KuIAQAAAAk1pgNAGZWKyl65md/l4F5QKefwBv3DsLRKwrdWd7IAAAAgOQaswFA0qt1+h6Ae4o7FD65/1Gk6TVmNuDJwGb2QknRaxHeOpxBAgAAAEkyJgOAmc2Q9DeRpmM6faIfdWNkuVnSewcp//HI8g5JPy17gAAAAEBCJSIAmNnzzOxfzGx5jL6rJf1C0rxI8+fc/UB//d39h5LuijR9zMxeVqL2pyVdEu3r7mHu5AQAAAAkQFJuBFYr6TpJ15nZI5J+Lmmt8sf1t0uaqPzx/n8g6TKdHlx+ojPvCVDsWkn3SpoiqUbSbWZ2s6TvSzogaaGkN0q6ILLNbZK+MZwvCgAAAEiapASAqDWFRxw3SXrHYJ/Su/tmM3ulpP9SPgRUSXpD4dGfn0t6nbtzhw0AAACMK4k4BEjSk5K+KWlLjL69yn86/yJ3v8bdjw22gSS5+z2SVkm6WdKJEt2ekvQ+SS929844dQEAAICxJBF7ANx9u6SrJMnMpkk6R9J0SVMlTZbUJemQpMclPTDU6/K7e5uk15tZs6SLJc2R1Chpj6RNku51dx/WFwMAAAAkWCICQJS771f+EJxKvsYR5Q8HAgAAAFIlKYcAAQAAABgBBAAAAAAgRQgAAAAAQIoQAAAAAIAUIQAAAAAAKUIAAAAAAFKEAAAAAACkCAEAAAAASBECAAAAAJAiBAAAAAAgRQgAAAAAQIoQAAAAAIAUIQAAAAAAKUIAAAAAAFKEAAAAAACkCAEAAAAASBECAAAAAJAiBAAAAAAgRQgAAAAAQIoQAAAAAIAUyY72AFBZJsk8TK1qZZTxMMWyygarVeNVyuZyQWpVW7WyPd1Baln1RKk3TK36KQs1ORfoP3LySi05cSJIqap5y1XrnUFqHWleqJqWniC19h6doIMds4dd5+iOLh3d2RRgRFJfbr9mr1gUpFb2xEQtXXBukFrH6xdJ08N8FpSZOFMNc8L8LjZMnamJ3b1BajXNWaRs7eEgtXLzlqurO8yfzp7aFh0/GuZr7M1NVI8WBqlV29OkyYtWB6k1acFyTa5vC1KratkKbfMw710TZ85WU2+Y9xtVVYWpA4wgAsA40lAz4Yy2eZkFwern6nuUsTBvdLm6blVlqsPUqpmsTKBafdmnZdmaMLUOP62q+jCTx6q2x9Uwb1mQWp3bH1NtXZg/op1bfi4/djhILR12HV6/MUipvkkXyg8/OfxCubnyrsPDryPJq5rVe3hvkFrK9ap926NBSlX31en4zjC1sq1SZu/aILVUtVzH234dpFS254Bq9z0WpJY0UbmdG8KUmpmRHd4cpJRNWKZMx/Ywtaqmq/OJXwWpVdX5lM5q2h+k1oH1/6OjG7cGqXV8+wT94Rf+OXb/jdv3KTfAh1brN20JMSxgxHAIEAAAwAAGmvwDYxEBAAAAAEgRAgAAAACQIgQAAAAAIEUIAAAAAECKEAAAAACAFCEAAAAAAClCAAAAAABShAAAAAAApAgBAAAAAEgRAgAAAACQIgQAAAAAIEUIAAAAAECKEAAAAACAFCEAAAAAAClCAAAAAABShAAAAAAApAgBAAAAAEgRAgAAAACQIgQAAAAAIEUIAAAAAECKEAAAAACAFCEAAAAAAClCAAAAAABShAAAAAAApAgBAAAAAEgRAgAAAACQIgQAAAAAIEUIAAAAAECKEAAAAACAFCEAAAAAAClCAAAAAABShAAAAAAApAgBAAAAAEgRAgAAAACQIgQAAAAAIEUIAAAAAECKEAAAAACAFCEAAAAAAClCAAAAAABShAAAAAAApAgBAAAAAEiR7GgPAOGsW7dutIcwJmTMtGLxvMq/kGXkoWrVNsu9N0gpm75AmcN7gtSqnTlDVU9vDFKrrvlstXZZkFodk+aqu7tv2HUm103X7oYwb5O1kxtV1TctSK2GqVM1oS7M92rC3KVqyHYHqVXfOk+9DUfC1GqZo67snCC1JrS2qKe9IUgtr66WWZjv/ZRZM1RTFeb32hvna1cmzLiydfWaOWNhkFp1U1qkY/slSd0veadUXTvwBu5Sie9vo6QL/ii/3HuiS7/5x08MWOr893xC2dq6ks+v37Z34LEA4xgBAKmTc1emrrHf56wmzCRBkmzizGAThUxzq6yqJkit6mP7g9U6vvEeWWtLkFqHNt2jqj07gtTKnsgpu2/DsOvkZj5Xjb4pwIgk92Wq6ggTljL1C3W87b4gtbJ9HfJd9wep1XfoUb32/e8MUst7u1U1/UNhasmVqb06TK2cZIEm2v/6xsuV9TAfE/TNu0zVx7YEqZU5UaNX/V2Y770kZSbkg2+sD6livmdma+t07U23DdjnsfXrY9UC0ohDgAAAAIAUIQAAAAAAKUIAAAAAAFKEAAAAAACkCAEAAAAASBECAAAAAJAiBAAAAAAgRQgAAAAAQIoQAAAAAIAUIQAAAAAAKUIAAAAAAFKEAAAAAACkCAEAAAAASBECAAAAAJAiBAAAAAAgRcoOAGbWZGYXm9n7zexmM9tsZjkz88LjzuEOysxWmNlnzOxBM9trZl1mts3M7jCzt5hZ4zBqn2dmXzCztWZ2wMw6zWyLmX3fzF5nZrXDHT8AAACQVNlyOpvZJklLJVklBmNmWUkfl/RnkqqKnp5feFwq6aNmdo27/6KM2hMkfV7Stf08vbjweKWkdWZ2lbs/XP5XAAAAACRbuXsAlqlCk/+CGyVdr1OTf5e0XtLdknZG+s2T9GMzuzROUTOrlnSbTp/890haK+mXkvZE2s+WdLeZrRnKFwAAAAAk2VDPAWhXflL+D5LeIOmh4Q7EzN4n6epI092SVrj7Kne/yN3nSXqxpN2F57OSvmNm82OU/5yk34+s3yppobuvcfcLJc2S9FpJRwvPN0q6fTiHGgEAAABJVG4AuFLSCknNhUn5+9z933Vq4jwkZjZV0sciTQ9JutTdN0f7uftPJV0oqaPQ1CTpk4PUXi7p7ZGm2yVd4e67InVz7v5tSS+V1FdoniPpg+V/NQAAAEBylRUA3P1b7r7J3T3wON4pqTmyfp27nygxhid0+qT/SjNbMEDtj+jUuQ49kt5Wavzu/r+Sbog0vdfMGgYZOwAAADBmJOUyoFdElu9z9/sH6f8VSV2F5Yyk1/TXqXDs/ysjTd+LfvJfwpciyxOV3ysAAAAAjAujHgDMbJGkVZGm2wfbxt0PSro30vSKEl0vkDS5zNrrJT0ZozYAAAAw5ox6AJD0zKL1X8XcLtrvGRWsXVwDAAAAGLOSEABWFq0/HnO7aL8mM5szSO0eSduGUHuZmRXfkwAAAAAYk5IQABZElvt06jKfg9k+QJ3+2na5e24ItWsltcbcDgAAAEi0JASApshyu7v3lex5uiNF6/1dsz9a+3AZY4pTGwAAABhzkhAAJkaWj5exXXHf/ibplawtSTKza83sATN74NChQ2W8BAAAADDykhAAqiPLvWVsV9y3up8+lawtSXL3G9z9PHc/b/LkyaW6AQAAAImQhABwLLJcV8Z2xX2P9dOnkrUBAACAMScJAaAjslzOXXeL+7aPcG0AAABgzElCANgXWZ5gZnFPuC2+Ms/+QWqXcyWf4r4HytgWAAAASKwkBICNRevzY24X7ZeTtHmQ2lPKCBfR2m3ufjjmdgAAAECiJSEAPFa0fm7M7aL9trl7f1f5Ka4d966+0drrY24DAAAAJF4SAsADOv0k24tibndhZPnOEn3uKloftLaZ1Uk6P0ZtAAAAYMwZ9QBQ+OT+R5Gm15jZgCfsmtkLJS2KNN1aovYOSQ9Gmq4yMxtkSJdLqh+sNgAAADAWjXoAKLgxstws6b2D9P94ZHmHpJ/GrL1U0utLdTSzWkl/Fmn6tbtzCBAAAADGjUQEAHf/oU4/XOdjZvay/vqa2aclXRLt6+7dA5T/iqQtkfUvmtmz+6mblXSDpJWR5o8MNnYAAABgLMmW09nMrpd0fT9P1USWLzSzrn76fMPd3zpA+Wsl3StpSqHebWZ2s6TvK38ZzoWS3ijpgsg2t0n6xkBjdvceM3uj8nsJaiVNlnS3md0o6SfKX+N/uaTrJK2JbPoldy8+hwAAAAAY08oKAIX+tYP0sRJ9qgfayN03m9krJf2X8iGgStIbCo/+/FzS69w9N8h45O73mNmVkr6u/E2+6iS9o/Doz7ckvWewugAAAMBYk4hDgE5y93skrZJ0s6QTJbo9Jel9kl7s7p1l1L5V0jMk/UBSX4luj0u6yt2vdPdSfQAAAIAxq6w9AO7+CUmfqMhITr1Gm6TXm1mzpIslzZHUKGmPpE2S7nV3H2LtxyVdZmYzlL+M6Gzl9wY8LWmdu/92+F8BAAAAkFzlHgI0Ytz9iPKHA1Wi9l5J361EbQAAACDJEnUIEAAAAIDKIgAAAAAAKUIAAAAAAFKEAAAAAACkSGJPAgYwBmSqpKqawfvFUN+6UFNzA94uJLZps89SY03vsOv0tS7W/rqBbjQeX+2cxcrt7u8eieVrmLtYzYcPBKlVNX2+drVvDlKrdmJTkDppMW3eArmGdFG7MzQsnaNGOx6kllUxNQDGOxviFTXRj1WrVvktt9wy2sMAYstkMlq5cuWQt8/1DX+S/bta7XtlZkFqdW99UNZb6lYi8e3esFmHntoVYERSX+1kNelIkFqqadD8NSuClHKZ6l5wpSRp46bNyuUGvbciKiBjphWL5gx5+64Hb5e62mP3f3LGc+UZJvpA0q1evfpBdz8vdF0OAQJSjMkeovh5GD25Ef4wjsk/kG4EAAAAACBFCAAAAABAihAAAAAAgBQhAAAAAAApQgAAAAAAUoQAAAAAAKQIAQAAAABIEQIAAAAAkCIEAAAAACBFCAAAAABAihAAAAAAgBQhAAAAAAApQgAAAAAAUoQAAAAAAKQIAQAAAABIEQIAAAAAkCIEAAAAACBFCAAAAABAihAAAAAAgBQhAAAAAAApQgAAAAAAUoQAAAAAAKQIAQAAAABIEQIAAAAAkCIEAAAAACBFCAAAAABAihAAAAAAgBQhAAAAAAApQgAAAAAAUoQAAAAAAKQIAQAAAABIEQIAAAAAkCIEAAAAACBFCAAAAABAihAAAAAAgBQhAAAAAAApQgAAAAAAUoQAAAAAAKQIAQAAAABIEQIAAAAAkCIEAAAAACBFCAAAAABAihAAAAAAgBQhAAAAAAApQgAAAAAAUoQAAAAAAKQIAQAAAABIEQIAAAAAkCIEAAAAACBFCAAAAABAihAAAAAAgBTJjvYAAIyudevWjfYQwqteIFVLcpfMBu9fqt+5y1V77skuLotRa6B+vYO9XtxxSXoiYC1t2Dj49qi49U/sHPrGk9acWo7zMxH35wbAuMQeAADjV9wJTox+cSb/sfsFHFfQWhgfQv7cABiXCAAAAABAihAAAAAAgBQhAAAAAAApQgAAAAAAUoQAAAAAAKQIAQAAAABIEQIAAAAAkCIEAAAAACBFCAAAAABAihAAAAAAgBQhAAAAAAApQgAAAAAAUoQAAAAAAKQIAQAAAABIEQIAAAAAkCIEAAAAACBFCAAAAABAihAAAAAAgBQhAAAAAAApQgAAAAAAUoQAAAAAAKRI2QHAzJrM7GIze7+Z3Wxmm80sZ2ZeeNw5hJoXR7Yv5/GjIbzWeWb2BTNba2YHzKzTzLaY2ffN7HVmVltuTQAAAGCsyJbT2cw2SVoqySoznMoxswmSPi/p2n6eXlx4vFLSOjO7yt0fHsHhAQAAACOirAAgaVlFRnGmO2L2eyBOJzOrlnSbpN+PNPdI2iDpiPJf18xC+9mS7jazC9z9kZjjAAAAAMaEcgPASe2SHpL0YOHxfknPDDUod39JqFoFn9Ppk/9bJb3b3XdJkpllJF0u6cuSmiQ1SrrdzFa6e3vgsQAAAACjptwAcKXyE/7N7u4nG83srUFHFZCZLZf09kjT7ZKuiI7f3XOSvm1mT0m6W1KVpDmSPijpYyM4XAAAAKCiyjoJ2N2/5e6bopPnMeAjOhV0eiS9rdT43f1/Jd0QaXqvmTVUeHwAAADAiBnXlwEtHPv/ykjT904e9jOAL0WWJ0p6afCBAQAAAKNkXAcASRdImhxZv32wDdx9vaQnI02vCD0oAAAAYLSM9wBQfGLyr2JuF+0X7ORmAAAAYLQlMgCY2U1mtsnM2s2sy8x2m9mvzOxvzOzcMkqtjCz3SNoWc7vHI8vLzKyqjNcEAAAAEiuRAUDS1cpfm3+ipFpJrZKeL+nDkh40sx+Z2dwYdRZElncVrvYTx/bI8snXBwAAAMa8od4HoNLalP+0/pjy1+VfXvj3pD+Q9LCZvcTd7x+gTnSbw2W8/pGi9cZSHc3sWp26u3DH6tWrN5XxOgAAAEAp8ytRNCkBwCXdKekmSXe4+9PRJws36rpA+Wvyn7yh1xTlb9b1LHd/qkTdiZHl42WMp7hvyQDg7jfo9EuHAgAAAImViEOA3P0ud/89d/9a8eS/8HzO3e+SdImkv488NUPSZwYoXR1Z7i1jSMV9q/vtBQAAAIwxiQgAcXneByX9LNJ8pZnNKLHJschyXRkvVdz3WL+9AAAAgDFmTAWAiOin/hlJl5bo1xFZLueOvsV928vYFgAAAEissRoAfqn8ZT1PWlGi377IcjlX8inue6CMbQEAAIDEGpMBwN17dPqkfFqJrhsjy1PMrOTJvEWiZ1y3ufvhMoYHAAAAJNaYDAAF0cN0Sl3h57Gi9bh39Y3ebGx97BEBAAAACTcmA0DhJmDRa/y3leh6V9H6RTFq10k6P9J0Z1mDAwAAABJsTAYASW8oWr+nv07uvkPSg5Gmq8zMBql9uaT6yPqt5Q8PAAAASKYxFwDMbLmkD0eadkn69QCb3BhZXirp9QPUrpX0Z5GmX7s7hwABAABg3Bj1AGBmrzazz5rZnBh9L1L+HgDNkeaPunvfAJt9RdKWyPoXzezZ/dTOKn9H35WR5o8MNiYAAABgLDF3j9/Z7HpJ1/fzVI2kk4fWuKTufvp8w93f2k/NayR9tbDdryXdLelR5S/heUz5Y/1XSnq5pIuLNv+au78xxrhfKOmnkmoLTV3K7xn4ifLX+F8u6TpJayKbfcnd3zVYbQAAAGAsyQ6hf+0gfaxEn+oY2z2v8BhMTtLnJP15jL5y93vM7EpJX1f+6kF1kt5RePTnW5LeE6c2AAAAMJaM+iFAkh6R9B1JO2P07ZL075Ke4+4fcvfeuC/i7rdKeoakH0gqdcjQ45KucvcrBzmsCAAAABiTyjoEqNLMbJakVcrf2Guq8sf6d0o6pPz1+B8q3ARsuK8zQ9KFkmYrvzfgaUnr3P23w60NAAAAJFmiAgAAAACAykrCIUAAAAAARggBAAjEzLaZmRce14z2eAAAwNgRmUO4mV1cydciAKBfZvaJoh9EN7NPlbF9XdG211RwuAAAICHMrNnMrjOzW8xsk5kdNLNeMztmZrvN7F4z+7qZvcfMnmtmVaM95rQhAKAc7zWz1tEeRKXxST4AAOUzsyoz+5Ck3ZL+RdIfS1omabKkKuUvxd4q6bmSrpL0D5LulbTfzM4flUFXwEh+kj9U5d4HAOnWIOnjkt422gMBAADJYWbVkr4t6VVFT3VL2ixpv/L3fJoqaalOv2fUJOVDAkYIewBQrjeb2bLRHgQAAEiUv9Lpk//HJf2JpEnuvtrdf8/dL3b31ZIaJZ0v6dPKhwOMMAIA4jgiaW9hOSvpM6M4FgAAkCBmNlPS+yNNayU9291vcffjxf3dvcfd73P36yWtkPRiSVtHZrSQCACIp0vSJyPrrxlPx+oBAIBh+UNJ1ZH1D7r7kTgbet5P3Z09ASOIAIC4/lXSE5H1z4Z+ATNrMrN3mtkPzGy7mXWa2VEze9zMvmVmfxL3SgH9nchrZlkz+2Mz+26hZnvh+a+Z2cUn+0uaHyn11X6uhnSyX9yv66WFKyFsMbPjZnbIzB4ys78xs5ayvkkAACTPiqL1X1XiRULNE8xsQdHf9AUxX7/kRULM7JoS84NflJhHbIv5mlkzu8LM/tvMnjSzLjPbb2a/NrPrzaw5Tp1inASMWNy9x8yul3RzoekiM3uZu/8gRH0zu0rS5yVN6+fpRklLJL1O0kfN7E3u/psy6y9SfuzPGe5Yy3jNqZK+qvwnI1F1kp5ReLzDzF4b6vsIAMAoaCpat9AvUOl5QhKZ2Xzl5y7PK3qqVvmTqc+X9C4zu8zdHyinNgEA5bhF0gckPauw/tdm9iN3zw2nqJn9ufInAkXtVf7EoGpJK5X/5VZh+edm9mp3vyPmS0yR9DNJCwrrbcrvzTDlr0QgSQclnax3kfKTdElaJ2lX3K8lokHSjyWdW/SaLmmVTl3tYKKk75nZs9x93RBeBwCA0ba/aP1SSf8ZqvgIzBNC2KVT84g/iLTfr/wco9ieQerNkHSTpHmF9R2Stin/9Z4jaUKk3x1mtsrd22KP1t158DjjIekTyk9WXVJbpP1FkXaXdHWJ7euK+l1Tot9Li/rtlPQKSZmiWu+UdCzS75CkuQOMf1uk79HCv2slXSzJIv0ykhYOsG2/447xmvsK/z4s6aKiflWS3iWpN9L/J6P9f86DBw8ePHgM5aH8nu7o3/Jdyp8EHKJ28HmC8h8IRmsuiDmWWPODotoXl/G1Rrc7OY+4U9Izi/rVSfpUUf8vl/N95RwAlMXdfybpJ5GmvzKz2lL9B2JmWeVvFHJSm6QL3P02j+xVcPcud/+S8pcX6ys0T1J+V2AcjcpP/l/o7nd64benUDvn7k8OZfyDmKZ86n+hu98VfcLd+9z9i8qHrJNeVNjVBwDAWHOH8pPjk2ZJ+o2Z/axwzP6aodztdwTnCUk0TdJ/S7rE3R+KPlH4eq+XdGOk+XVmVh+3OAEAQ/Fh5dOmlD9h9h1DrPNHOrVrS5Le7e7bSnV2959I+ufo9mY2r1T/Ite6+9HyhzhkOUlXuXvHAH3+SfkrLEn5w5FeUPFRAQAQmLt3S7pa0olIs0n6fUlfVH5veHvhxNUvmNnlMU9eHcl5QtK0K7+XoXeAPn8bWZ4g6ZlxixMAULZCEr0l0vTnQzwL/Y8iyzskfSfGNp/TqfBRpfxuwME84iN/MtDP3H3TQB0KgeThSNOqio4IAIAKcfdfKv9B1mMlutQrf9Lqnyr/977NzL45yM1FR2qekET/4e79nTvwO56/dGr0uP/Y8wgCAIbqLyT1FJanKr9XoFzRs9r/J3poTinuvl3SIyVqlHJ3uQMLIO4l0J6KLE+qwDgAABgR7v6g8ieoXq784StdA3Svk3SlpMfM7N0l+ozUPCGJKjqPIABgSNx9q/L3Bjjp3WbWGnf7wrGACyJNj5To2p+1keUlMfqPxt0F456Jfyyy3FCJgQAAMFIK59bd6u6vUH5CeoHydwn+hvJX7SmWlfSPZvauaOMIzxOSqKLzCAIAhuOvlD9GTcr/0H2ijG0nF63vK2PbaN/iOv0ZyWP/TzoxeJczBL9uMgAAo8XdT7j7Pe7+eXe/2t2XK3/u4EeVv0pP1N+Z2ZzI+kjOE5KoovMIAgCGzN33KX+s3UlvMrPlMTcvvnJQdxkvHf2lqCvZ65Rh3acAAACE4e473P1TklZLip4rVyvpzUXrUZWcJ6QOAQDD9TmduplFVtJnYm53uGi9sb9OJUTvOFhcBwAAJJy775L0tqLmCyLLh4ueS8I8oexLmSYVAQDDUrjM5ScjTa82s/NjbHdMUmekaXEZLxvtu7eM7QAAQHLcJSl6uexZJxcqPE8o3ptQHbPuWD2c6AwEAIRwg6QtkfW/LdWxyIOR5efH2aBwU5DnRJoeiPla5YgeMsRx+QAAVEDhqj7RAFB8zftKzROKzw2cEqPuQuWvtR9H9GpFiZxHEAAwbO7eI+n6SNOFki6LsWn0DrmXmFlLjG0uU/6yo/3VCCV6Rn3su+oBAID4zGyKpBmRpt1FXSoyTygcvRC9xv45Meq+OkafkxI/jyAAIJRv6/SUHedcgBt16tP2akmfHaizmdVK+utI03ZJPy5jjHE9HVleWoH6AACMG2Z2UeET8nL9qU6fi/606PlKzhOic5bXDlJ3mqQPDNSnSOLnEQQABFHYjRe9GdhAd/Y7uc02STdFmq42s4+Z2Rm7y8xsgvJ3Hz4r0vxJd+8b2ogHFN3l+Nqiy5IBAIDTvVjSZjO7qRAGBjzsxcyqzOwDyl8O9KR2Sf8e7VfhecJ3I8svMrOrSox1hqTbJcXZ+3BSdB7xJjNL3LkD2dEeAMYPd/+5mf1Y0qVlbPYeSRdJWlRY/0tJl5nZ15S/PFi1pHMlXavTbwjyfXe/cZhDLuVbyocZk9Qq6Qkz+63y1xX+3fkB7v6qCr0+AABjTVbS1YXHTjO7S9J9knYof7hNlaSZyv9Nf43OPKn3ve7+tM70HlVmnvAt5Q9fnldYv8nMLpX0PeX/3k9V/qpEb5bUrPzeiRWS4nwo+E1Jf1JYPkfSU4V5xEGdOj9gr7tfG6NWRRAAENqHlf8kINZJL+5+1MwukvQjSasKzc/R6SfwFPtPSa8fziAHGdOjZvZx5W90Jkk1kp5bqdcDAGCMK/6Ufa6kNxQeg+mQ9G53/7f+nqzUPMHdjxU+9f+h8jcztQHG/FtJfyzpoYFqRmr/j5l9RdJbCk0Nkl5Y1G17nFqVwiFACMrdH5Z0c5nbPCXp2ZL+QgPf6W+zpCslXe7uQ7lDXjlj+qTyyf/fJD2m/BUDuKEYAABn+kvlLwDy98pPluMcnrtL0uclrSg1+T+pUvMEd79b+b0L95fo0q78/Y6e7+7Fdy4erPZblT8h+ebCuDp0+tWBRpXlD90GksHMMsr/kq+SNF35S4LtlXS/u28czbEBAIDBFY7HXylpifJX+Zmo/N/zduVPkH1U0hM+hElopeYJZna28nv7pxfGuV3Szwv3Ixh3CAAAAABAinAIEAAAAJAiBAAAAAAgRQgAAAAAQIoQAAAAAIAUIQAAAAAAKUIAAAAAAFKEAAAAAACkCAEAAAAASBECAAAAAJAiBAAAAAAgRQgAAAAAQIr8f24TD0/6upTFAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "width_scale = 35.8/45.2\n", "fig, ax = plt.subplots(1, 1, figsize=(15*width_scale, 10*150/250))\n", "imin, imax = 234-2, 301+2\n", "gridding = 2\n", "sisal = 30\n", "idsal = 31.2\n", "\n", "vmax = 0.01 / transports[sals < sisal].sum() / m2lat * 24 * 365 * 4\n", "print ('oo', transports.sum(), transports[sals < sisal].sum(), transports[sals >= sisal].sum())\n", "\n", "c, xedge, yedge, im = ax.hist2d(lons[sals < sisal], \n", " depths[sals < sisal],\n", " weights=transports[sals < sisal],\n", " bins=[np.array(f_lats[imin:imax:gridding]), new_ws], cmap='Greens_alpha', \n", " vmax=vmax, density=True)\n", "print (c.max()/vmax)\n", "\n", "v2max = 0.01 / otransports[osals <= idsal].sum() / m2lat * 24 * 365 * 4\n", "print ('oo', otransports.sum(), otransports[osals < idsal].sum(), otransports[osals >= idsal].sum())\n", "c, xedge, yedge, im2 = ax.hist2d(olons[osals < idsal], \n", " odepths[osals < idsal],\n", " weights=otransports[osals < idsal],\n", " bins=[np.array(f_lats[imin:imax:gridding]), new_ws], cmap='Oranges_alpha', \n", " vmax=v2max, density=True,)\n", "print (c.max()/v2max)\n", "\n", "\n", "\n", "vmax = 0.01 / transports[sals >= sisal].sum() / m2lat * 24 * 365 * 4\n", "c, xedge, yedge, im = ax.hist2d(lons[sals >= sisal], \n", " depths[sals >= sisal],\n", " weights=transports[sals >= sisal],\n", " bins=[np.array(f_lats[imin:imax:gridding]), new_ws], cmap='Blues_alpha', \n", " vmax=vmax, density=True,)\n", "\n", "print (c.max()/vmax)\n", "\n", "v2max = 0.01 / otransports[osals > idsal].sum() / m2lat * 24 * 365 * 4\n", "c, xedge, yedge, im2 = ax.hist2d(olons[osals >= idsal], \n", " odepths[osals >= idsal],\n", " weights=otransports[osals >= idsal],\n", " bins=[np.array(f_lats[imin:imax:gridding]), new_ws], cmap='Purples_alpha', \n", " vmax=v2max, density=True)\n", "\n", "print (c.max()/v2max)\n", "\n", "\n", "\n", "ax.invert_yaxis()\n", "ax.set_ylim(150, 0)\n", "ax.invert_xaxis()\n", "ax.text(48.37, 170, 'North')\n", "ax.text(48.15, 170, 'South')\n", "ax.set_xticks([])\n", "pc = draw_patches_vic(fmask2, imin, imax, np.array(v_lats), np.array(w_depths[0]))\n", "ax.add_collection(pc);\n", "\n", "print ((mymesh.e2f[0, imin:imax, jgrid].sum()).values)\n", "fig.savefig('victoria_final.png')" ] }, { "cell_type": "code", "execution_count": 96, "id": "42e84310-627d-4b91-8c61-d6b9af3ee450", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[26.59728865 26.28684983 25.70331479 24.64073198 22.81323256 19.95870258\n", " 16.10792044]\n" ] } ], "source": [ "print (e3t[30:23:-1].values)" ] }, { "cell_type": "code", "execution_count": 99, "id": "636f726c-b7ed-43dd-be28-c8361c7abedf", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array(19.79981227)" ] }, "execution_count": 99, "metadata": {}, "output_type": "execute_result" } ], "source": [ "(e3t[23]+e3t[22]).values" ] }, { "cell_type": "code", "execution_count": 105, "id": "27c6634a-847a-4215-a968-1e76ff773aa1", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array(14.79823468)" ] }, "execution_count": 105, "metadata": {}, "output_type": "execute_result" } ], "source": [ "(e3t[21:15:-1].sum()).values" ] }, { "cell_type": "code", "execution_count": 108, "id": "8a6599f1-784e-4116-a80c-857e7c477e2c", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array(15.09450132)" ] }, "execution_count": 108, "metadata": {}, "output_type": "execute_result" } ], "source": [ "(e3t[:15].sum()).values" ] }, { "cell_type": "code", "execution_count": 109, "id": "d8aeeb8b-d115-4cf5-8e61-c70a83af938e", "metadata": {}, "outputs": [], "source": [ "def sumthemup(array):\n", " newarray = np.empty((10))\n", " newarray[0] = array[0:15].sum()\n", " newarray[1] = array[15:21].sum()\n", " newarray[2] = array[22] + array[23]\n", " newarray[3:10] = array[24:31]\n", " return newarray" ] }, { "cell_type": "code", "execution_count": 18, "id": "5adbc405-b8d0-4dd2-8103-3da8648470ae", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([ 0. , 15.09625559, 31.10103509, 50.96323669,\n", " 67.05207424, 86.96747037, 109.73706598, 134.34593444,\n", " 160.02956159, 186.30527815, 212.89655735])" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "new_ws = np.empty((11))\n", "new_ws[0] = w_depths[0, 0]\n", "new_ws[1] = w_depths[0, 15]\n", "new_ws[2] = w_depths[0, 22]\n", "new_ws[3:] = w_depths[0, 24:32]\n", "new_ws" ] }, { "cell_type": "code", "execution_count": null, "id": "5311adfa-a60b-4285-8ec8-f22323894703", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "py39", "language": "python", "name": "py39" }, "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.2" } }, "nbformat": 4, "nbformat_minor": 5 }