Skip to content

Instantly share code, notes, and snippets.

@willirath

willirath/parcels_NEMO_tracer_vert-axis.ipynb

Last active August 28, 2023 16:04
Show Gist options
  • Save willirath/88689100f3ce573d50c551cf9c4e116f to your computer and use it in GitHub Desktop.
Save willirath/88689100f3ce573d50c551cf9c4e116f to your computer and use it in GitHub Desktop.
Download ZIP
Parcels C-Grid tracer vert axis supplementary
Display the source blob
Display the rendered blob
Raw
parcels_NEMO_tracer_vert-axis.ipynb
{
"cells": [
{
"cell_type": "markdown",
"id": "1329ad66-cbaa-47a7-8cd2-9b8bca6322f5",
"metadata": {},
"source": [
"# Reproducible example"
]
},
{
"cell_type": "markdown",
"id": "ccc0f4bd-31f1-4536-b51f-1d5315dbe215",
"metadata": {},
"source": [
"Download example data including tracer file from https://doi.org/10.5281/zenodo.3634490"
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "9a8af5b2-46d6-417f-b3c8-6966025a86ad",
"metadata": {},
"outputs": [],
"source": [
"# !curl -O https://zenodo.org/record/3634491/files/NEMO_GYRE_test_data_all_files.v2020.02.03.1.zip"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "2958d4c2-03d1-4de5-b602-be91c4d1428b",
"metadata": {},
"outputs": [],
"source": [
"# !mkdir -p GYRE_5d\n",
"# !unzip NEMO_GYRE_test_data_all_files.v2020.02.03.1.zip -d GYRE_5d"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "67de1594-3eec-47c0-840d-45ea42d883a7",
"metadata": {},
"outputs": [],
"source": [
"import xarray as xr\n",
"\n",
"from pathlib import Path"
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "99a12220-ca02-438a-8b54-f33743a1db47",
"metadata": {},
"outputs": [],
"source": [
"data_files = sorted(Path(\"GYRE_5d/\").glob(\"*.nc\"))\n",
"\n",
"for f in data_files:\n",
" try:\n",
" ds = xr.load_dataset(f, decode_cf=False)\n",
" ds.time_counter.attrs[\"calendar\"] = \"360_day\"\n",
" ds.to_netcdf(f)\n",
" except Exception as e:\n",
" print(e)"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "d9420cb6-2e13-4277-a7bf-d0c6a988a38b",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAjMAAAGwCAYAAABcnuQpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA8Z0lEQVR4nO3deXBUdbr/8U+HJQZIghGzaQRUwlXjjsOiNaAsIyKiaF0dKYWqGa6joqLj6FXvHeMCuS4s13EFLYxeEatUlBlXEAI6yMhaN26og2guW5QRAjgmLOf3h0X/OlunO5zT/Zxz3q+qrppe8uT5dHfidzrfh2/EcRxHAAAAPpWR7gYAAAAOBYsZAADgayxmAACAr7GYAQAAvsZiBgAA+BqLGQAA4GssZgAAgK91THcDXjtw4IA2b96s7OxsRSKRdLcDAAAS4DiOdu3apeLiYmVkxP/sJfCLmc2bN6ukpCTdbQAAgHaoqanR0UcfHfcxgV/MZGdnS/r5ycjJyUlzNwAAIBF1dXUqKSmJ/nc8nsAvZg7+aSknJ4fFDAAAPpPIFhE2AAMAAF9jMQMAAHyNxQwAAPA1FjMAAMDXWMwAAABfYzEDAAB8jcUMAADwNRYzAADA11jMAAAAX2MxAwAAfI3FDAB4bO/+A3GvAzg0gT+bCQDSyXEcLV3/nR5b8pXWb9ulvgXZuv7c4zX0hPyEzpwB0DY+mQEAj+zdf0DvfVaric+v0tqaHfqxYb/W1uzQxOdX6b3PavmEBnAJixkA8EinDhl6bMlXcpzGtzuO9HjVV+rUgV/BgBv4SQIAD63ftqvl27e2fDuA5LGYAQAP9S3Ibvn2wpZvB5C8tC5mysvLFYlEGl0KCwuj9zuOo/LychUXFysrK0tDhgzRJ598ksaOAaC51qaV9u4/oOvPPV5N9/lGItJ1Q45v9HVMPAHtl/ZpppNOOkmLFi2KXu/QoUP0fz/44IOaPn26nn32WZWWlur+++/X8OHDtX79emVn8/9qAKRfvGmlTh0yNPSEfM2+qp8er/pK67fuUt/CbF03pPE0ExNPwKGJOE7TrWmpU15ertdee03r1q1rdp/jOCouLtbkyZN1++23S5Lq6+tVUFCgBx54QNdcc02LNevr61VfXx+9XldXp5KSEu3cuVM5OTme5AAQTnv3H9DS9d9p4vOrGm3yjUSk2Vf10+C+R6pThwzt3X+g0Wbf2OuJ1gDCpq6uTrm5uQn99zvtPyFffvmliouL1bt3b11xxRXasGGDJOnrr7/W1q1bNWLEiOhjMzMzNXjwYC1fvrzVehUVFcrNzY1eSkpKPM8AIJwSnVZquhiJvc7EE3Do0vpT0r9/fz333HN65513NHv2bG3dulWDBg3S9u3btXXrVklSQUFBo68pKCiI3teSO+64Qzt37oxeampqPM0AINzcmFZi4gk4NGndMzNy5Mjo/z755JM1cOBAHXfccaqsrNSAAQMkqdnfix3Hifs35MzMTGVmZnrTMAA00bcgW2trdjS/PYlpJTdqAGFm6vPLrl276uSTT9aXX34ZnWpq+ilMbW1ts09rACAdkplW8rIGEHamFjP19fX67LPPVFRUpN69e6uwsFALFy6M3t/Q0KClS5dq0KBBaewSQJjEG5mOnVY645ju6tq5g844prtmX9UvOs3kVo226gBhltY/M916660aPXq0jjnmGNXW1ur+++9XXV2dxo8fr0gkosmTJ2vq1Knq06eP+vTpo6lTp6pLly668sor09k2gJBIZGQ6EolocN8jNezE//+J8d79B5Iau26rRqJ1gLBK62j2FVdcoWXLlun777/XkUceqQEDBui+++7TiSeeKOnnH9577rlHTz31lH744Qf1799fjz32mMrKyhL+HsmMdgHAQW6MTLs1ds34NsIomf9+p3UxkwosZgC01yWP/bXFjblnHNNdr153dspquFkH8Atf/TszAGCVpbFrxreB1rGYAYBWuHFIpFsHTXJgJdA6FjMAQivedJAbh0S6ddBkonWYdkJYsWcGQCg5jqP3PquNOx108DFtHRIZr44bNRLpJZEagJ+wATgGixkATSUzHeTGIZFuHTTZWh2mnRBEbAAGgDiSOdzRjUMi3TposrU6HFaJsOMdDiCULE0ZWakB+BWLGQChZGnKyEoNwK9YzAAIHbcOd7Ry0CSHVSLs2AAMILDibbxNZMqorRqJ1nGjRqryAFYwzRSDxQwQTomMKie6yGhr3DmRRcah1Eh1HsACFjMxWMwA4WPlkEhLB00yvg2/YTQbQKi5MapspYa1XgCLePcCCCQr485BGwEHLGIxAyCQrIw7B20EHLCIxQwAX/L6YMagHTTpVi+ARWwABuA7qTqYMWgHTbrVC5AKTDPFYDEDBEuqD2YM2kGTbvUCeI1pJgCBleqDGYN20KRbvQCW8M4E4DuWJnus9GIpD5BqLGYA+I6lyR4rvVjKA6QaixkAvmLpYEYrvVjKA6QDG4ABmOOngxmDdtCkG3kANzDNFIPFDOAvfjyYMWgHTbqRBzhULGZisJgB/CNoBzOSB2g/RrMB+FLQDmYkD5AavPMAmGJpxNhKDUu9MLoNi1jMADDF0oixlRqWemF0GxaxmAGQcmE5mJE8LfcCuI0NwABSKmwHM5KHiSe0D9NMMVjMAHaE9WBG8jDxhOQxzQTApLAezEgeJp7gLd5VAFLKykQNU0b2ewESxWIGQEpZmahhysh+L0CiWMwASBkrhypaOpjRSg1rvQDJYAMwAFdxMCN53OoF4cY0UwwWM0DqcDAjedzMg3Dz5WKmoqJCd955p2666SbNnDlTkjRhwgRVVlY2elz//v21YsWKhOuymAFSw9JBhlZ6IY+3dRBsvhvNXrlypWbNmqVTTjml2X3nn3++tmzZEr28+eabaegQQFssHWRopRfyeFsHOCjt75jdu3dr3Lhxmj17tg4//PBm92dmZqqwsDB6ycvLi1uvvr5edXV1jS4AUsPSWK+VXsjjbR1AMrCYuf766zVq1CgNGzasxfurqqqUn5+v0tJSTZw4UbW1tXHrVVRUKDc3N3opKSnxom0ALbA01mulF/J4WweQ0ryYmTdvntasWaOKiooW7x85cqReeOEFLV68WNOmTdPKlSt13nnnqb6+vtWad9xxh3bu3Bm91NTUeNU+EDp+OcgwaAczkqflmsBBadsAXFNTo379+undd9/VqaeeKkkaMmSITjvttOgG4Ka2bNminj17at68eRo7dmxC34cNwIA7/HaQYdAOZiQP005h44tpptdee02XXHKJOnToEL1t//79ikQiysjIUH19faP7DurTp49++9vf6vbbb0/o+7CYAQ6dXw8yDNrBjORh2ilMfDHNNHToUFVXV2vdunXRS79+/TRu3DitW7euxYXM9u3bVVNTo6KiojR0DISXXw8yDNrBjORh2gkt65iub5ydna2ysrJGt3Xt2lVHHHGEysrKtHv3bpWXl+vSSy9VUVGRNm7cqDvvvFM9evTQJZdckqaugfCyNMVipYalXsiDMDO7vO3QoYOqq6s1ZswYlZaWavz48SotLdWHH36o7Gx2uwOpZmmKxUoNS72QB2FmajFTVVUV3fyblZWld955R7W1tWpoaNA333yjZ599llFrIA2CdpAheWz34lYehIeZ4wy8wgZgIHFhOsiQPOHIA//yxTRTqrCYARITxoMMyROOPPAnFjMxWMwAbXNjFNZKDUu9kMd+L7DLF6PZAOywcgihpYMMrdSw1EvQ8iA4eLUBSLIzTsuIse1egpYHwcBiBoAkO+O0jBjb7iVoeRAMLGaAkOAgQ/KQp+Xr8D82AAMhwEGG5CEPE09+wzRTDBYzCLtUH/wXb5w20Tpu1CAPeZh48jemmQBEcZAhecjTdh34G68kEAKWJkes9EIe271YygP7WMwAIWBpcsRKL+Sx3YulPLCPxQwQAEyxkIc83uSBP7ABGPA5pljIQx5v8yA9mGaKwWIGQcYUC3nIk5o8SD2mmYCQYIqFPORJTR7YxisF+JylqQ8rNSz1Qh7bvTDtFAwsZgCfszT1YaWGpV7IY7sXpp2CgcUM4GPJTH14XcdKDUu9kMd2L27lQfqxARjwgXgbIROZ+mirRqJ13KhBHvIEMQ/cxzRTDBYz8LtERkcT/aXf1vhpIr/0D6UGecgT5DxwF4uZGCxm4GdujI5aqWGpF/LY7iVoedA+jGYDAeHG6KiVGpZ6IY/tXoKWB97jVQCMszJ+ykiu7V7IY78XeIfFDGCclfFTRnJt90Ie+73AOyxmgDRrOv4Zez3R0VE3asSr40YN8pAnrHngPTYAA2nEwX/kIU848iB5TDPFYDEDq5KZkmhtdNSNGsn04kYN8pAnrHmQHKaZAB9IZkqi6S/Cg9fdqJFMHTdqkIc8Yc0D7/AMA2lkadLCSi/ksd0Lebytg/ZhMQOkkaVJCyu9kMd2L+Txtg7ah8UMkCbJTEl4WcNSL+Sx3Qt5vK2D9mMDMOCheJsGJVsH5QXt4D/ykMdvedAY00wxWMwgXfx4UF4iv/TJQx7yeJcH/x+LmRgsZpAObo1qulHHSg1LvZDHdi/kgcRoNpB2bo1qulHHSg1LvZDHdi/kQbJ4BgGPWBr5tFLDUi/ksd0LeZAMM4uZiooKRSIRTZ48OXqb4zgqLy9XcXGxsrKyNGTIEH3yySfpaxJIgqWRTys1LPVCHtu9kAfJMLGYWblypWbNmqVTTjml0e0PPvigpk+frkcffVQrV65UYWGhhg8frl27WMnChqYjlwevJzOq2VqNZOq4UYM85CGP7TxoXdo3AO/evVtnnHGGHn/8cd1///067bTTNHPmTDmOo+LiYk2ePFm33367JKm+vl4FBQV64IEHdM0117RYr76+XvX19dHrdXV1KikpYQMwXBe2g/LIQx7ypD9PmPhqmmn8+PHKy8vTjBkzNGTIkOhiZsOGDTruuOO0Zs0anX766dHHjxkzRt27d1dlZWWL9crLy3XPPfc0u53FDNyU6HRCvFHNZCYcWqvjRg3ykIc8/soTFr6ZZpo3b57WrFmjioqKZvdt3bpVklRQUNDo9oKCguh9Lbnjjju0c+fO6KWmpsbdpgGF96A88pCHPOnPg+Y6pusb19TU6KabbtK7776rww47rNXHNf1ozXGcuB+3ZWZmKjMz07U+gdZYmXBg6sN2L+Sx3UvQ8oRV2pZ6q1evVm1trc4880x17NhRHTt21NKlS/XII4+oY8eO0U9kmn4KU1tb2+zTGiAdrEw4MPVhuxfy2O4laHnCKm2LmaFDh6q6ulrr1q2LXvr166dx48Zp3bp1OvbYY1VYWKiFCxdGv6ahoUFLly7VoEGD0tU2ICm56QTrNSz1Qh7bvZDHfi9hlfYNwLFiNwBL0gMPPKCKigrNmTNHffr00dSpU1VVVaX169crOzuxlSrHGaC94m3WkzgojzzkIY+9PEHiq2mmWE0XM47j6J577tFTTz2lH374Qf3799djjz2msrKyhGuymEF7cFAeechDHr/mCQrfLma8wGIGyXJjRNKtMUsrvZDHdi/ksd2LpTx+4pvRbMAiN0Yk3RqztNILeWz3Qh7bvVjKE1ThTQ7EYWnM0kov5LHdC3ls92IpTxCxmAFaYGnM0kov5LHdC3ls92IpTxCxmEEoNR1zjL2ezIhka3XcqJFMHfKQhzzkSSRPULEBGKHjt4PlgnZQHnnIQ5705vELpplisJhBrGSmAeKNSCZax40a8XohD3nIQ5725PEDppmAViQzDdD0h77pfRyURx7ykMeveYImuMmAVliaKrBSw1Iv5LHdC3ls9xLWaScWMwgdS1MFVmpY6oU8tnshj+1ewjrtxGIGoZLMVIHXdazUsNQLeWz3Qh7bvbiVx4/YAIxAirfRLmgHy5GHPOQhj9t5LGCaKQaLmfAJ48Fy5CEPecjjdp50YzETg8VMuLgxmmilhqVeyGO7F/LY7iVoeVKF0WyElhujiVZqWOqFPLZ7IY/tXoKWxyJ/dg3EYWW8kZFP272Qx3Yv5LHfiyUsZhA4VsYbGfm03Qt5bPdCHvu9WMJiBr7TdLww9nqio4lu1IhXx40a5CEPechjMY9FbACGr3CwHHnIQx7ypD9PKjDNFIPFTHAkswu/tdFEN2ok04sbNchDHvKQx2IerzHNhEBKZhd+0x+0g9fdqJFMHTdqkIc85CGPxTyW2OsIiMPSTn4rvZDHdi/ksd0LebytkyosZuArlnbyW+mFPLZ7IY/tXsjjbZ1UYTED30hmF76XNSz1Qh7bvZDHdi/k8bZOKrEBGKbE25Qm2TqILWgHy5GHPOQhj5t5DhXTTDFYzPiHHw9iS+SXCnnIQx7yhDXPoWAxE4PFjD/s3R+sg9jIY7sX8tjuhTy2e3ErT1sYzYbvuDUK6EYdKzUs9UIe272Qx3Yv5PEeixmYYWmk0EoNS72Qx3Yv5LHdC3m8xWIGZlgaKbRSw1Iv5LHdC3ls90Ieb7GYQUo1Hek7eD2ZUcDWaiRTx40a5CEPechDnvg1UoUNwEiZsB3ERh7ykIc8Yc5zqJhmisFixoZEd7/HGwVMZgd9a3XcqEEe8pCHPOSJn8cNTDPBnER3vzf9QYi9nswO+tbquFGDPOQhD3nIE79GqrGYQcpY2UHPVIHtXshjuxfy2O7FUp5UYjGDlLGyg56pAtu9kMd2L+Sx3YulPKnEYgYp4cbudys1LPVCHtu9kMd2L+Txtk4qsQEYrmlrMxgHsZGHPOQhT3DyeM0300xPPPGEnnjiCW3cuFGSdNJJJ+mPf/yjRo4cKUmaMGGCKisrG31N//79tWLFioS/B4uZ1OAgNvKQhzzkCU+eVPDNYubPf/6zOnTooOOPP16SVFlZqYceekhr167VSSedpAkTJmjbtm2aM2dO9Gs6d+6svLy8hL8Hixnv7d1v5+AyK72Qx3Yv5LHdC3ns95IKvhnNHj16tC644AKVlpaqtLRUU6ZMUbdu3Rp98pKZmanCwsLopa2FTH19verq6hpd4K1kxgG9rGGpF/LY7oU8tnshj/1erDHT9f79+zVv3jzt2bNHAwcOjN5eVVWl/Px8lZaWauLEiaqtrY1bp6KiQrm5udFLSUmJ161DwRsptFLDUi/ksd0LeWz3ErQ81qR9MVNdXa1u3bopMzNTv/vd7zR//nydeOKJkqSRI0fqhRde0OLFizVt2jStXLlS5513nurr61utd8cdd2jnzp3RS01NTaqihFrQRgqt1LDUC3ls90Ie270ELY81aV/M9O3bV+vWrdOKFSt07bXXavz48fr0008lSZdffrlGjRqlsrIyjR49Wm+99Za++OILvfHGG63Wy8zMVE5OTqMLDl3TUbzY68mM8bVWx40aydQhD3nIQx7yxO/FTxLeAHzLLbfovvvuU9euXbVs2TINGjRIHTt2dL2hYcOG6bjjjtNTTz3V4v19+vTRb3/7W91+++0J1WMD8KGzdHBZ0A5iIw95yEMei3ks8GSaqVOnTvq///s/FRQUqEOHDtqyZYvy8/NdaTjW0KFDVVJSomeffbbZfdu3b9dRRx2lWbNm6eqrr06oHouZQ5PMzvd4Y3yJ1nGjRrxeyEMe8pCHPPHzWOHJNFOvXr30yCOPaOnSpXIcRx9++KGWLVvW4iVRd955p95//31t3LhR1dXVuuuuu1RVVaVx48Zp9+7duvXWW/Xhhx9q48aNqqqq0ujRo9WjRw9dcsklCX8PHJpkdr43/UFoeh8HsZGHPOQhj/08fpTw34keeugh/e53v1NFRYUikUirC4pIJKL9+/cnVHPbtm266qqrtGXLFuXm5uqUU07R22+/reHDh+uf//ynqqur9dxzz2nHjh0qKirSueeeq5deeknZ2f7dpORHQduFTx7bvZDHdi/ksd1LECeVEpHwUuziiy/W1q1bVVdXJ8dxtH79ev3www/NLv/4xz8S/ubPPPOMNm7cqPr6etXW1mrRokUaPny4JCkrK0vvvPOOamtr1dDQoG+++UbPPvsso9ZpELRd+OSx3Qt5bPdCHtu9BHFSKRFJf67UrVs3LVmyRL17927077nEXuA/QduFTx7ykIc85Gm9RtAkvAG4rq4uugGnrX9V19JGWzYAty1ou/DJQx7ykIc8tieVEuHJNFPsBFNGRkaLT4rjOEntmUkFFjPxBW0XPnnIQx7ykMf+pFIiPJlmWrx4cfRcpCVLlmjx4sXNLgdvh38EbRc+echDHvKQJ1iTSolIeJpp8ODBLf5v+J+VHfRMFdjuhTy2eyGP7V4s5Qmidi/XduzYoXfffVf/8z//o+eee67RBf5iZQc9UwW2eyGP7V7IY7sXS3mCqF2LmT//+c865phjNHLkSE2aNEk33XRT9DJ58mSXW4SX3Nj9bqWGpV7IY7sX8tjuhTze1gmihDcAxyotLdUFF1ygqVOnqkuXLl705Ro2AMffUCYltoPejRpt1XGjBnnIQx7yhD1PUHgyzRSra9euqq6u1rHHHtvuJlMl7IuZRMYBpcR+CA+lRqJ13KhBHvKQhzxhzRMkni9mxo4dqyuuuEL/+q//2u4mUyXMi5m9+xMfB/SyhqVeyGO7F/LY7oU89nsJEk9GsxcsWBC9jBo1Sn/4wx9UXl6uV155pdF9CxYsOOQAcEcy44Be1rDUC3ls90Ie272Qx34vYZXwaPbFF1/c7LZ777232W3W/tG8sAvaSKGVGpZ6IY/tXshju5eg5QmrhJd6Bw4cSOjCQsaWoI0UWqlhqRfy2O6FPLZ7CVqesOJzK59rOooXez2ZMb7W6rhRI5k65CEPechDnvi9oLl2bQCWpPfee08zZszQZ599pkgkon/5l3/R5MmTNWzYMLd7PCRB3gBs6eCyoB3ERh7ykIc8FvOEiefTTI8++qhuvvlmXXbZZRo4cKAkacWKFXr55Zc1ffp0TZo0qX2deyCoi5lkdr7HG+NLtI4bNeL1Qh7ykIc85ImfJ2w8mWaKVVFRoRkzZujFF1/UjTfeqBtvvFFz587VjBkzNHXq1HY1jeQks/O96Q9C0/s4iI085CEPeeznQeva9SzV1dXp/PPPb3b7iBEjVFdXd8hNITFB24VPHtu9kMd2L+Sx3QuTSt5q12Lmoosu0vz585vd/vrrr2v06NGH3BQSE7Rd+OSx3Qt5bPdCHtu9MKnkrXYtZk444QRNmTJFo0aN0v3336/7779fF154oaZMmaKTTjpJjzzySPQCb7i1892NOlZqWOqFPLZ7IY/tXsiDZLVrA3Dv3r0TKx6JaMOGDUk35Sa/bwCOtxmMg9jIQx7ykCdYefD/eT7N5Cd+XsxYObgsaAexkYc85CGPxTxozPPFzC233NJysUhEhx12mPr06aOLLrpIeXl5yZZ2nV8XM3v32zi4zI0alnohj+1eyGO7F/J4WweNeT6avXbtWj3zzDOaNWuWli5dqqqqKs2ePVvPPPOM3nvvPd188806/vjj9cknn7QrAOwcXOZGDUu9kMd2L+Sx3Qt5vK2D9mvXMzxmzBgNGzZMmzdv1urVq7VmzRpt2rRJw4cP169//Wtt2rRJv/zlL1v9BAeJsTIOyIik7V7IY7sX8tjuxVIetF+7FjMPPfSQ7rvvvkYf++Tk5Ki8vFwPPvigunTpoj/+8Y9avXq1a42GkZVxQEYkbfdCHtu9kMd2L5byoP3atZjZuXOnamtrm93+3XffRf/RvO7du6uhoeHQugu4puN4sdcTHeVzo0a8Om7UIA95yEMe8sBL7doAPG7cOH344YeaNm2azjrrLEUiEX300Ue69dZbNWjQID3//POaN2+eHn74Ya1atcqLvhNmdQOwpYPLgnYQG3nIQx7yWMuD5Hk+zbR7927dfPPNeu6557Rv3z5JUseOHTV+/HjNmDFDXbt21bp16yRJp512WtIB3GRxMZPMzvfWRvncqJFML27UIA95yEOeMOZB+3g+zdStWzfNnj1b27dv19q1a7VmzRpt375ds2bNUteuXSX9vIhJ90LGqmR2vjf9QTh43Y0aydRxowZ5yEMe8oQxD7zX8VC+uFu3bjrllFPc6iVUgrYL30oNS72Qx3Yv5LHdS9DywFssG9MkaLvwrdSw1At5bPdCHtu9BC0PvMViJg3c2Pnu1u55K72Qx3Yv5LHdC3ns9wJvcTaTR9raCMZBbOQhD3nIE6w8cBcHTcZIx2LG0sFlQTuIjTzkIQ95LOaB+1jMxEj1YmbvfjsHl1mpYakX8tjuhTy2eyEPUsnz0Wy3PPHEEzrllFOUk5OjnJwcDRw4UG+99Vb0fsdxVF5eruLiYmVlZWnIkCHmD69MZozP6zpWaljqhTy2eyGP7V7IA6vS+kodffTR+q//+i+tWrVKq1at0nnnnacxY8ZEFywPPvigpk+frkcffVQrV65UYWGhhg8frl27bI/DBW2kkDy2eyGP7V7IY7sXxq6DIa2LmdGjR+uCCy5QaWmpSktLNWXKFHXr1k0rVqyQ4ziaOXOm7rrrLo0dO1ZlZWWqrKzUjz/+qLlz56az7TYFbaSQPLZ7IY/tXshjuxfGroPBzGdo+/fv17x587Rnzx4NHDhQX3/9tbZu3aoRI0ZEH5OZmanBgwdr+fLlrdapr69XXV1do4sXmo7jHbyezBhfazWSqeNGDfKQhzzkIU/8GrAt7RuAq6urNXDgQP3000/q1q2b5s6dqwsuuEDLly/X2WefrU2bNqm4uDj6+H/7t3/TN998o3feeafFeuXl5brnnnua3e7mBmArB5cF7SA28pCHPOSxmAfp4atppoaGBn377bfasWOHXnnlFT399NNaunSpduzYobPPPlubN29WUVFR9PETJ05UTU2N3n777Rbr1dfXq76+Pnq9rq5OJSUlri1mEt39Hm+ML5kd9K3VcaMGechDHvKQJ34epI9vppkkqXPnzjr++OPVr18/VVRU6NRTT9V///d/q7CwUJK0devWRo+vra1VQUFBq/UyMzOj01EHL25KdPd70x+E2OvJ7KBvrY4bNchDHvKQhzwcEBkE5l4xx3FUX1+v3r17q7CwUAsXLoze19DQoKVLl2rQoEFp7NDODnqmCmz3Qh7bvZDHdi+W8sC+tC5m7rzzTr3//vvauHGjqqurddddd6mqqkrjxo1TJBLR5MmTNXXqVM2fP18ff/yxJkyYoC5duujKK69MZ9tmdtAzVWC7F/LY7oU8tnuxlAf2pXUxs23bNl111VXq27evhg4dqr/97W96++23NXz4cEnSbbfdpsmTJ+u6665Tv379tGnTJr377rvKzk7fm9CN3e9WaljqhTy2eyGP7V7I420d2Jf2DcBe8+I4Aw5iIw95yEOe4OSBTb6aZvKaV2czJfJDyEFs5CEPechjOw/sYjETg4Mm01/DUi/ksd0LeWz3ErQ8sM1Xo9lBk8xIodd1rNSw1At5bPdCHtu9BC0PgoNX2wOWRgqt1LDUC3ls90Ie270ELQ+CgcWMByyNFFqpYakX8tjuhTy2ewlaHgQDixmXJTMK2HQssOl9fjuIjTzkIQ95LOZB8LEB2ANhPIiNPOQhD3ks5oF/Mc0UIx2LGSn+KGAyu/Bbq+NGjWR6IQ95yEMev+WBvzHNZEDTH6bY68nswm+tjhs1kqlDHvKQhzx+y4Pw4FVPE0s7+a30Qh7bvZDHdi/kQZixmEkTSzv5rfRCHtu9kMd2L+RBmLGYSQM3duG7tZPfSi/ksd0LeWz3Qh6EHRuA08RvB7G5UYM85CEPedzOg+BimimG1cWMlPgvhHQfxBa0g+XIQx7yBCMPgo3FTAzLi5l4khlNtF7DUi/ksd0LeWz3YikPgo/R7ABIZjTReg1LvZDHdi/ksd2LpTxALN4xhlkZb2Tk03Yv5LHdC3m8rQNILGZMszLeyMin7V7IY7sX8nhbB5BYzJiV6Ghi0xHF2OvJjDe2VseNGuQhD3nI095egESwAdgwDpYjD3nIE9Y8ANNMMfy8mJFaH01MZhog3nhjonXcqEEe8pCHPBwSiUQxzRQgTX+wD15PZhqgtRrJ1HGjBnnIQx7yJNoLkAzeOT5maarASg1LvZDHdi/ksd8LkCgWMz5maarASg1LvZDHdi/ksd8LkCgWMz7l1jSAG3Ws1LDUC3ls90Ie+70AyWADsI+F9WA58pCHPMHIA8TDNFOMIC9mpMR/qQTlYDnykIc8wcgDtIXFTIygL2biSWZE0ssalnohj+1eyGO7F7fyAIlgNBuSgnewnJUalnohj+1eyAOkBu+8gLM0ZmmlF/LY7oU8tnth7BoWsZgJOEtjllZ6IY/tXshjuxfGrmERi5kAS2ZEsum45MHrbtRIpo4bNchDHvKkPw+QSmwADjg/HSwXtIPyyEOesOYB3MA0U4ywL2ak+COSiU4nuFEjXi9u1CAPeciT/jyAW5hmQiNNf8nEXk90OsGNGvHquFGDPOQhT/rzAOnAOxBmJhyY+rDdC3ls92IpD5BqLGZgZsKBqQ/bvZDHdi+W8gCpltbFTEVFhc466yxlZ2crPz9fF198sdavX9/oMRMmTFAkEml0GTBgQJo6Dh4/TkkEbeqDPOTxWx7AmrRuAD7//PN1xRVX6KyzztK+fft01113qbq6Wp9++qm6du0q6efFzLZt2zRnzpzo13Xu3Fl5eXkJfQ82ALfNT1MSQZv6IA95/JYHSBXfTjN99913ys/P19KlS/XLX/5S0s+LmR07dui1115rV00WM4nxw5RE0KY+yEMev+UBUsm300w7d+6UpGafulRVVSk/P1+lpaWaOHGiamtrW61RX1+vurq6Rhe0zQ9TEkGb+iAPefyWB7DKzDvUcRzdcsstOuecc1RWVha9feTIkXrhhRe0ePFiTZs2TStXrtR5552n+vr6FutUVFQoNzc3eikpKUlVhMCyNCVhpYalXshju5eg5QEsMrOYmTRpkv73f/9XL774YqPbL7/8co0aNUplZWUaPXq03nrrLX3xxRd64403Wqxzxx13aOfOndFLTU1NKtoPNEtTElZqWOqFPLZ7CVoewCITi5kbbrhBCxYs0JIlS3T00UfHfWxRUZF69uypL7/8ssX7MzMzlZOT0+iC9nNrusGNOlZqWOqFPLZ7CVoewKq0bgB2HEc33HCD5s+fr6qqKvXp06fNr9m+fbuOOuoozZo1S1dffXWbj2cD8KFLdLqhrU2DidRxo0ZbdchDHvIwqQT7fDPNdN1112nu3Ll6/fXX1bdv3+jtubm5ysrK0u7du1VeXq5LL71URUVF2rhxo+688059++23+uyzz5Sd3fZHoyxm3JHoL8l4I59t1XGjRqJ1yEMe8rR8HbDCN4uZ1v6fwJw5czRhwgT985//1MUXX6y1a9dqx44dKioq0rnnnqv77rsv4Y29LGa8l8zIp5c1LPVCHtu9kAewzzej2Y7jtHiZMGGCJCkrK0vvvPOOamtr1dDQoG+++UbPPvssE0rGJDPy6WUNS72Qx3Yv5AGChXc4XGFpbNRKL+Sx3Qt5gOBgMQNXWBobtdILeWz3Qh4gOFjM4JAlM/LZdPzz4HU3aiRTx40a5CGPH/MAQWTqbCYvsAE4Nfx0UF7QDv4jD3kSzQP4iW+mmVKBxUzqxBv5THTawo0a8XpxowZ5yOPHPIDf+GaaCcHS9Jdm7PVEpy3cqBGvjhs1yEMeP+YBgox3OlLGysQGUyy2eyGPt3WAIGIxg5SxMrHBFIvtXsjjbR0giFjMICXcmLawUsNSL+Sx3YulPECQsQEYKcPBf+Qhj7d5gCBhmikGixlbEvmlz8F/5CEPh0QCLGZisJjxh2RGWL2uY6WGpV7IY78XIGgYzYbvJDPC6nUdKzUs9UIe+70AYcZPCcywNMJqpYalXshjvxcgrFjMwAxLI6xWaljqhTz2ewHCisUMTEhm9LTpGGrT+/x28B95yMPYNXBo2AAMM8J48B95yJNoHiBsmGaKwWLGX+KNniYz9dFaHTdqJNMLecjjRh4gjJhmgm+5cdhevDpu1EimDnnI40YeAPHx0wJfsTQ5YqUX8tjuhUklwHssZuArliZHrPRCHtu9MKkEeI/FDHzDjakPtyZHrPRCHtu9MKkEpAYbgOErfjv4z40a5AlHHgCNMc0Ug8VM8CT6H6B0H/wXtIMMyeNdHgDNsZiJwWImXJIZhbVew1Iv5PG2DoDmGM1GaCUzCmu9hqVeyONtHQCHhp80BI6VcVpGjG33YikPgEPDYgaBY2WclhFj271YygPg0LCYQaBwkCF5LOYB4C02ACNwOMiQPJbyAGgfpplisJgJJw4yJI+FPADaj2kmhB4HGZLHQh4AqcFPHELH0hSLlRqWeglaHgDeYzGD0LE0xWKlhqVegpYHgPdYzCBUOMjQdi9BywMgNdgAjNAJ60GG5PEuDwD3Mc0Ug8UMWhK2gwzJ410eAN7wzWKmoqJCr776qj7//HNlZWVp0KBBeuCBB9S3b9/oYxzH0T333KNZs2bphx9+UP/+/fXYY4/ppJNOSuh7sJhBsoJ2kKGVGpZ64YBIwD7fjGYvXbpU119/vVasWKGFCxdq3759GjFihPbs2RN9zIMPPqjp06fr0Ucf1cqVK1VYWKjhw4dr1y6mCeCNoB1kaKWGpV44IBIIlo7p/OZvv/12o+tz5sxRfn6+Vq9erV/+8pdyHEczZ87UXXfdpbFjx0qSKisrVVBQoLlz5+qaa65JR9sIAUtjvVZ6IQ8Aq0z934+dO3dKkvLy8iRJX3/9tbZu3aoRI0ZEH5OZmanBgwdr+fLlLdaor69XXV1dowuQLEtjvVZ6IQ8Aq8wsZhzH0S233KJzzjlHZWVlkqStW7dKkgoKCho9tqCgIHpfUxUVFcrNzY1eSkpKvG0cgePHgwyDdjCjpTwA7DMzzXT99dfrjTfe0AcffKCjjz5akrR8+XKdffbZ2rx5s4qKiqKPnThxompqapr9mUr6+ZOZ+vr66PW6ujqVlJSwARhJ8dNBhkE7mNFSHgDp45tppoNuuOEGvfbaa1q2bJl69+4dvX3Dhg067rjjtGbNGp1++unR28eMGaPu3bursrKyzdpMM6G9/HCQYdAOZrSUB0B6+WaayXEcTZo0Sa+++qoWL17caCEjSb1791ZhYaEWLlwYva2hoUFLly7VoEGDUt0uQsYPBxkG7WBGS3kA+Edap5muv/56zZ07V6+//rqys7Oj+2Byc3OVlZWlSCSiyZMna+rUqerTp4/69OmjqVOnqkuXLrryyivT2TpgZqKGKSNv6wCwL63/N+SJJ57Qzp07NWTIEBUVFUUvL730UvQxt912myZPnqzrrrtO/fr106ZNm/Tuu+8qO5uJA6SXlYkapoy8rQPAvrT/mamly4QJE6KPiUQiKi8v15YtW/TTTz9p6dKl0WknIF2sHGQYtIMZLeUB4B8mNgB7iQ3A8AoHM/o/DwC7fDfN5CUWM/ASBzP6Nw8A21jMxGAxg3TgYEb7vQCwzTej2UBQcTCj/V4ABAc/9YBHLI0YW6lhrRcAwcBiBvCIpRFjKzWs9QIgGFjMAB4I88GMlvIACAc2AAMeCePBjJbyAPA3pplisJhBOoXpYEZLeQD4H9NMgBFhOpjRUh4A4cJPP5BGliZ7rPTCpBKAZLGYAdLI0mSPlV6YVAKQLBYzQJoE7WBGKzUAhA8bgIE08tvBjJYOmgQQbEwzxWAxA+v8cjCjpYMmAQQfi5kYLGbgZ1YOZrR00CSAcGA0GwgIKwczWjpoEgCa4jcHYJyVcWdLI+AAEIvFDGCclXFnSyPgABCLxQxgWFgPmgSAZLABGDAubAdNAoDENFMjLGYQBGE5aBIADmKaCQiYsBw0CQDtwW8QwOcsTRkxqQQgHVjMAD5nacqISSUA6cBiBvCxVE8ZMakEwCI2AAM+l6opIyaVAKQS00wxWMwgDLyeMmJSCUCqMc0EhIzXU0ZMKgGwjN8yQAgwqQQgyFjMACHApBKAIGMxAwScG1NGTCoBsIzFDBBwnTpkaOgJ+Zp9VT+dcUx3de3cQWcc012zr+qnoSfkN9usG+vg9WRqAECqMc0EhERbU0aJjF4zqQQgVZhmAtBMvCmjvfsP6L3PajXx+VVaW7NDPzbs19qaHZr4/Cq991lto09o4tUEgHTgNxGApEavAcAafkMBkMToNQD/SutiZtmyZRo9erSKi4sViUT02muvNbp/woQJikQijS4DBgxIT7NAwDF6DcCv0rqY2bNnj0499VQ9+uijrT7m/PPP15YtW6KXN998M4UdAuHA6DUAP+uYzm8+cuRIjRw5Mu5jMjMzVVhYmKKOgHCKHb3mkEgAfpPWxUwiqqqqlJ+fr+7du2vw4MGaMmWK8vPzW318fX296uvro9fr6upS0Sbge5FIRIP7HqlhJxZEb9u7/wALGQDmmd4APHLkSL3wwgtavHixpk2bppUrV+q8885rtFhpqqKiQrm5udFLSUlJCjsG/I3RawB+ZOYfzYtEIpo/f74uvvjiVh+zZcsW9ezZU/PmzdPYsWNbfExLn8yUlJTwj+YBAOAjyfyjeeb/zBSrqKhIPXv21JdfftnqYzIzM5WZmZnCrgAAQDr56jPk7du3q6amRkVFReluBQAAGJHWT2Z2796tr776Knr966+/1rp165SXl6e8vDyVl5fr0ksvVVFRkTZu3Kg777xTPXr00CWXXJLGrgEAgCVpXcysWrVK5557bvT6LbfcIkkaP368nnjiCVVXV+u5557Tjh07VFRUpHPPPVcvvfSSsrP5R7wAAMDPzGwA9gqnZgMA4D+cmg0AAEKDxQwAAPA1FjMAAMDXWMwAAABf89U/mtceB/c3c0YTAAD+cfC/24nMKQV+MbNr1y5J4owmAAB8aNeuXcrNzY37mMCPZh84cECbN29WdnZ2yk//PXguVE1NTWjHwsP+HJA/3PklngPyhzu/1P7nwHEc7dq1S8XFxcrIiL8rJvCfzGRkZOjoo49Oaw85OTmhfRMfFPbngPzhzi/xHJA/3Pml9j0HbX0icxAbgAEAgK+xmAEAAL7GYsZDmZmZuvvuu5WZmZnuVtIm7M8B+cOdX+I5IH+480upeQ4CvwEYAAAEG5/MAAAAX2MxAwAAfI3FDAAA8DUWMwAAwNdYzHjkiy++0JgxY9SjRw/l5OTo7LPP1pIlSxo9JhKJNLs8+eSTaerYXYnk//bbbzV69Gh17dpVPXr00I033qiGhoY0deyuqqqqFl/fSCSilStXRh8X1PdAovmD/B6QpDfeeEP9+/dXVlaWevToobFjxza6P6ivf6y2noMgvwd69erV7PX993//90aPCfJ7IJH8br3+gf8XgNNl1KhRKi0t1eLFi5WVlaWZM2fqwgsv1N///ncVFhZGHzdnzhydf/750euJ/muH1rWVf//+/Ro1apSOPPJIffDBB9q+fbvGjx8vx3H0pz/9Kd3tH7JBgwZpy5YtjW77z//8Ty1atEj9+vVrdHsQ3wOJ5A/6e+CVV17RxIkTNXXqVJ133nlyHEfV1dXNHhfE1/+gtp6DoL8HJOnee+/VxIkTo9e7devW7DFBfg/Ey+/q6+/Add99950jyVm2bFn0trq6OkeSs2jRouhtkpz58+enoUNvJZL/zTffdDIyMpxNmzZFH/Piiy86mZmZzs6dO1Pes9caGhqc/Px859577210e1DfA021lD/I74G9e/c6Rx11lPP000/HfVyQX/9EnoMgvwccx3F69uzpzJgxI+5jgvweaCu/m68/f2bywBFHHKETTjhBzz33nPbs2aN9+/bpqaeeUkFBgc4888xGj500aZJ69Oihs846S08++aQOHDiQpq7dk0j+Dz/8UGVlZSouLo5+3a9+9SvV19dr9erV6WrdMwsWLND333+vCRMmNLsviO+BplrKH+T3wJo1a7Rp0yZlZGTo9NNPV1FRkUaOHKlPPvmk2WOD+von8hwE+T1w0AMPPKAjjjhCp512mqZMmdLin1CC+h6Q4ud38/Xnz0weiEQiWrhwocaMGaPs7GxlZGSooKBAb7/9trp37x593H333aehQ4cqKytL7733nn7/+9/r+++/13/8x3+kr3kXJJJ/69atKigoaPR1hx9+uDp37qytW7emoWtvPfPMM/rVr36lkpKSRrcH9T3QVEv5g/we2LBhgySpvLxc06dPV69evTRt2jQNHjxYX3zxhfLy8iQF+/VP5DkI8ntAkm666SadccYZOvzww/XRRx/pjjvu0Ndff62nn346+pggvwfayu/q69/ej4/C6O6773Ykxb2sXLnSOXDggHPRRRc5I0eOdD744ANn9erVzrXXXuscddRRzubNm1ut//DDDzs5OTkpTJQcN/NPnDjRGTFiRLPv0alTJ+fFF19MdbSEJfocxKqpqXEyMjKcl19+uc36QXkPxGotvx/fA4nmf+GFFxxJzlNPPRX92p9++snp0aOH8+STT7Za3/rr7zjuPgdBfg+05OWXX3YkOd9//32r9a2/B9zM7+brzyczSZg0aZKuuOKKuI/p1auXFi9erL/85S/64YcfosedP/7441q4cKEqKyub7eY+aMCAAaqrq9O2bduarVYtcDN/YWGh/va3vzX62h9++EF79+41mf2gRJ+DWHPmzNERRxyhiy66qM36QXkPxGotvx/fA4nm37VrlyTpxBNPjN6emZmpY489Vt9++22rX2v99ZfcfQ6C/B5oyYABAyRJX331lY444ohWH2P5PeBmfjdffxYzSejRo4d69OjR5uN+/PFHSVJGRuMtSRkZGXH/Frp27Voddthhjf4UZYmb+QcOHKgpU6Zoy5YtKioqkiS9++67yszMbLavyJJEn4ODHMfRnDlzdPXVV6tTp05tPj4o74GD4uX343sg0fxnnnmmMjMztX79ep1zzjmSpL1792rjxo3q2bNnq19n/fWX3H0OgvweaMnatWslKZq1tcdYfg+4md/V1z/Zj5jQtu+++8454ogjnLFjxzrr1q1z1q9f79x6661Op06dnHXr1jmO4zgLFixwZs2a5VRXVztfffWVM3v2bCcnJ8e58cYb09z9oUsk/759+5yysjJn6NChzpo1a5xFixY5Rx99tDNp0qQ0d++uRYsWOZKcTz/9tNl9QX4PHBQvf9DfAzfddJNz1FFHOe+8847z+eefO7/5zW+c/Px85x//+IfjOOF4/dt6DoL8Hli+fLkzffp0Z+3atc6GDRucl156ySkuLnYuuuii6GOC/B5IJL+brz+LGY+sXLnSGTFihJOXl+dkZ2c7AwYMcN58883o/W+99ZZz2mmnOd26dXO6dOnilJWVOTNnznT27t2bxq7d01Z+x3Gcb775xhk1apSTlZXl5OXlOZMmTXJ++umnNHXsjV//+tfOoEGDWrwv6O8Bx4mf33GC/R5oaGhwfv/73zv5+flOdna2M2zYMOfjjz+O3h+G17+t58BxgvseWL16tdO/f38nNzfXOeyww5y+ffs6d999t7Nnz57oY4L8Hkgkv+O49/pHHMdx2vV5EQAAgAH8OzMAAMDXWMwAAABfYzEDAAB8jcUMAADwNRYzAADA11jMAAAAX2MxAwAAfI3FDAAA8DUWMwBM6tWrl2bOnOn59/nxxx916aWXKicnR5FIRDt27PD8ewJwFwdNAgi1yspKvf/++1q+fLl69Oih3NzcdLcEIEksZgCE2t///nedcMIJKisrS3crANqJPzMBSItdu3Zp3Lhx6tq1q4qKijRjxgwNGTJEkydPbvHx06dP18knn6yuXbuqpKRE1113nXbv3h29/9lnn1X37t31l7/8RX379lWXLl102WWXac+ePaqsrFSvXr10+OGH64YbbtD+/fslSUOGDNG0adO0bNkyRSIRDRkyJAXJAbiNxQyAtLjlllv017/+VQsWLNDChQv1/vvva82aNa0+PiMjQ4888og+/vhjVVZWavHixbrtttsaPebHH3/UI488onnz5untt99WVVWVxo4dqzfffFNvvvmmnn/+ec2aNUsvv/yyJOnVV1/VxIkTNXDgQG3ZskWvvvqqp5kBeIM/MwFIuV27dqmyslJz587V0KFDJUlz5sxRcXFxq18T+4lN7969dd999+naa6/V448/Hr197969euKJJ3TcccdJki677DI9//zz2rZtm7p166YTTzxR5557rpYsWaLLL79ceXl56tKlizp37qzCwkJvwgLwHIsZACm3YcMG7d27V7/4xS+it+Xm5qpv376tfs2SJUs0depUffrpp6qrq9O+ffv0008/ac+ePerataskqUuXLtGFjCQVFBSoV69e6tatW6PbamtrPUgFIF34MxOAlHMcR5IUiURavL2pb775RhdccIHKysr0yiuvaPXq1Xrsscck/fxpzEGdOnVq9HWRSKTF2w4cOHDIGQDYwWIGQModd9xx6tSpkz766KPobXV1dfryyy9bfPyqVau0b98+TZs2TQMGDFBpaak2b96cqnYBGMefmQCkXHZ2tsaPH68//OEPysvLU35+vu6++25lZGQ0+7RG+nnxs2/fPv3pT3/S6NGj9de//lVPPvlkGjoHYBGfzABIi+nTp2vgwIG68MILNWzYMJ199tk64YQTdNhhhzV77Gmnnabp06frgQceUFlZmV544QVVVFSkoWsAFkWc1v5IDQAptGfPHh111FGaNm2afvOb36S7HQA+wp+ZAKTF2rVr9fnnn+sXv/iFdu7cqXvvvVeSNGbMmDR3BsBvWMwASJuHH35Y69evV+fOnXXmmWfq/fffV48ePdLdFgCf4c9MAADA19gADAAAfI3FDAAA8DUWMwAAwNdYzAAAAF9jMQMAAHyNxQwAAPA1FjMAAMDXWMwAAABf+3/4k3rtvsbnBwAAAABJRU5ErkJggg==",
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"ds_meshmask = xr.open_dataset(\"GYRE_5d/mesh_mask.nc\").squeeze()\n",
"ds_meshmask.plot.scatter(x=\"glamf\", y=\"gphif\");"
]
},
{
"cell_type": "markdown",
"id": "23a4db7b-5c63-4a19-9636-89178c90bd99",
"metadata": {},
"source": [
"## Parcels Experiment\n",
"\n",
"(largely follows https://docs.oceanparcels.org/en/latest/examples/tutorial_nemo_3D.html )"
]
},
{
"cell_type": "code",
"execution_count": 6,
"id": "8f0ea17a-7fe4-4cc0-bf54-6a3bf297a5ef",
"metadata": {},
"outputs": [],
"source": [
"from datetime import timedelta as delta\n",
"from glob import glob\n",
"\n",
"import matplotlib.pyplot as plt\n",
"import xarray as xr\n",
"\n",
"from parcels import (\n",
" AdvectionRK4_3D,\n",
" FieldSet,\n",
" Field,\n",
" JITParticle,\n",
" ScipyParticle,\n",
" ParticleSet,\n",
" XarrayDecodedFilter,\n",
" logger,\n",
" Variable,\n",
")\n",
"\n",
"from pathlib import Path"
]
},
{
"cell_type": "code",
"execution_count": 7,
"id": "75003e45-188d-4a7c-88a3-c4a6d054fbdc",
"metadata": {},
"outputs": [],
"source": [
"# Add a filter for the xarray decoding warning\n",
"logger.addFilter(XarrayDecodedFilter())"
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "2e1785c5-4206-4bdf-9e37-4bc4a57464fd",
"metadata": {},
"outputs": [],
"source": [
"data_path = Path(\"GYRE_5d/\")\n",
"\n",
"tfiles = [data_path / \"GYRE_5d_00010101_00011230_grid_T.nc\", ]\n",
"ufiles = [data_path / \"GYRE_5d_00010101_00011230_grid_U.nc\", ]\n",
"vfiles = [data_path / \"GYRE_5d_00010101_00011230_grid_V.nc\", ]\n",
"wfiles = [data_path / \"GYRE_5d_00010101_00011230_grid_W.nc\", ]\n",
"mesh_mask = [data_path / \"mesh_mask.nc\", ]\n",
"\n",
"filenames = {\n",
" \"U\": {\"lon\": mesh_mask, \"lat\": mesh_mask, \"depth\": wfiles[0], \"data\": ufiles},\n",
" \"V\": {\"lon\": mesh_mask, \"lat\": mesh_mask, \"depth\": wfiles[0], \"data\": vfiles},\n",
" \"W\": {\"lon\": mesh_mask, \"lat\": mesh_mask, \"depth\": wfiles[0], \"data\": wfiles},\n",
"}\n",
"\n",
"variables = {\n",
" \"U\": \"vozocrtx\", \n",
" \"V\": \"vomecrty\", \n",
" \"W\": \"vovecrtz\",\n",
"}\n",
"\n",
"# Note that all variables need the same dimensions in a C-Grid\n",
"c_grid_dimensions = {\n",
" \"lon\": \"glamf\",\n",
" \"lat\": \"gphif\",\n",
" \"depth\": \"depthw\",\n",
" \"time\": \"time_counter\",\n",
"}\n",
"dimensions = {\n",
" \"U\": c_grid_dimensions,\n",
" \"V\": c_grid_dimensions,\n",
" \"W\": c_grid_dimensions,\n",
"}\n",
"\n",
"fieldset = FieldSet.from_nemo(filenames, variables, dimensions)"
]
},
{
"cell_type": "code",
"execution_count": 9,
"id": "e895db78-210b-4b16-a862-817ee14dc651",
"metadata": {},
"outputs": [],
"source": [
"temperature_field = Field.from_netcdf(\n",
" filenames=tfiles,\n",
" variable=(\"temperature\", \"votemper\"),\n",
" dimensions={\"lon\": \"nav_lon\", \"lat\": \"nav_lat\", \"depth\": \"deptht\", \"time\": \"time_counter\"},\n",
" interp_method=\"nearest\",\n",
")\n",
"\n",
"fieldset.add_field(temperature_field, \"temperature\")"
]
},
{
"cell_type": "code",
"execution_count": 10,
"id": "c384b51e-b6a3-45dd-8949-f967d6fcbbca",
"metadata": {},
"outputs": [],
"source": [
"class TemperatureParticle(JITParticle):\n",
" temperature = Variable(\"temperature\", initial=0)\n",
"\n",
"\n",
"def SampleTemperature(particle, fieldset, time):\n",
" particle.temperature = fieldset.temperature[particle]"
]
},
{
"cell_type": "code",
"execution_count": 11,
"id": "de860322-eb6b-4d92-bcd3-81de5ad6ea83",
"metadata": {},
"outputs": [],
"source": [
"pset = ParticleSet.from_line(\n",
" fieldset=fieldset,\n",
" pclass=TemperatureParticle,\n",
" size=1_000,\n",
" start=(-80.0, 30.0),\n",
" finish=(-55.0, 35.0),\n",
" # depth=float(fieldset.temperature.grid.depth[0] + 0.1), # this is working\n",
" depth=float(fieldset.temperature.grid.depth[0]), # this is working\n",
" # depth=float(ds_meshmask.gdept_1d.min().data[()] - 0.1), # this is NOT working\n",
")"
]
},
{
"cell_type": "code",
"execution_count": 12,
"id": "f7ba9185-e6f1-4a85-9492-f5c835c13ff4",
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"INFO: Compiled ArrayTemperatureParticleSampleTemperatureAdvectionRK4_3D ==> /tmp/parcels-1000/lib376f0b7129be8602e59bd4614ff599e6_0.so\n",
"WARNING: dt or runtime are zero, or endtime is equal to Particle.time. The kernels will be executed once, without incrementing time\n"
]
}
],
"source": [
"kernels = pset.Kernel(SampleTemperature) + pset.Kernel(AdvectionRK4_3D)\n",
"pset.execute(kernels, runtime=delta(days=0), dt=delta(hours=0))"
]
},
{
"cell_type": "markdown",
"id": "110a1215-8c51-4cad-9103-f23a85c365b9",
"metadata": {},
"source": [
"## Compare with direct selection"
]
},
{
"cell_type": "code",
"execution_count": 13,
"id": "7aba0b9b-7ae6-4323-80de-b4ea7fe39b15",
"metadata": {},
"outputs": [],
"source": [
"import pandas as pd"
]
},
{
"cell_type": "code",
"execution_count": 14,
"id": "e6bea8f6-df1b-4c9d-94ac-2ea4670b60bb",
"metadata": {},
"outputs": [],
"source": [
"def particle_to_pandas(particle):\n",
" return pd.DataFrame(\n",
" {\n",
" \"lon\": [particle.lon, ],\n",
" \"lat\": [particle.lat, ],\n",
" \"time\": [particle.time, ],\n",
" \"depth\": [particle.depth, ],\n",
" \"temperature\": [particle.temperature, ],\n",
" },\n",
" index=[particle.id, ],\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": 15,
"id": "a57da102-899f-4d88-86cc-849ad7c6bb97",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>lon</th>\n",
" <th>lat</th>\n",
" <th>time</th>\n",
" <th>depth</th>\n",
" <th>temperature</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>-80.000000</td>\n",
" <td>30.000000</td>\n",
" <td>0.0</td>\n",
" <td>4.975266</td>\n",
" <td>22.839523</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>-79.974975</td>\n",
" <td>30.005005</td>\n",
" <td>0.0</td>\n",
" <td>4.975266</td>\n",
" <td>22.839523</td>\n",
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <td>-79.949950</td>\n",
" <td>30.010010</td>\n",
" <td>0.0</td>\n",
" <td>4.975266</td>\n",
" <td>22.839523</td>\n",
" </tr>\n",
" <tr>\n",
" <th>3</th>\n",
" <td>-79.924925</td>\n",
" <td>30.015015</td>\n",
" <td>0.0</td>\n",
" <td>4.975266</td>\n",
" <td>22.869986</td>\n",
" </tr>\n",
" <tr>\n",
" <th>4</th>\n",
" <td>-79.899900</td>\n",
" <td>30.020020</td>\n",
" <td>0.0</td>\n",
" <td>4.975266</td>\n",
" <td>22.869986</td>\n",
" </tr>\n",
" <tr>\n",
" <th>...</th>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" </tr>\n",
" <tr>\n",
" <th>995</th>\n",
" <td>-55.100100</td>\n",
" <td>34.979980</td>\n",
" <td>0.0</td>\n",
" <td>4.975266</td>\n",
" <td>22.643122</td>\n",
" </tr>\n",
" <tr>\n",
" <th>996</th>\n",
" <td>-55.075075</td>\n",
" <td>34.984985</td>\n",
" <td>0.0</td>\n",
" <td>4.975266</td>\n",
" <td>22.643122</td>\n",
" </tr>\n",
" <tr>\n",
" <th>997</th>\n",
" <td>-55.050050</td>\n",
" <td>34.989990</td>\n",
" <td>0.0</td>\n",
" <td>4.975266</td>\n",
" <td>22.643122</td>\n",
" </tr>\n",
" <tr>\n",
" <th>998</th>\n",
" <td>-55.025025</td>\n",
" <td>34.994995</td>\n",
" <td>0.0</td>\n",
" <td>4.975266</td>\n",
" <td>22.643122</td>\n",
" </tr>\n",
" <tr>\n",
" <th>999</th>\n",
" <td>-55.000000</td>\n",
" <td>35.000000</td>\n",
" <td>0.0</td>\n",
" <td>4.975266</td>\n",
" <td>22.643122</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"<p>1000 rows × 5 columns</p>\n",
"</div>"
],
"text/plain": [
" lon lat time depth temperature\n",
"0 -80.000000 30.000000 0.0 4.975266 22.839523\n",
"1 -79.974975 30.005005 0.0 4.975266 22.839523\n",
"2 -79.949950 30.010010 0.0 4.975266 22.839523\n",
"3 -79.924925 30.015015 0.0 4.975266 22.869986\n",
"4 -79.899900 30.020020 0.0 4.975266 22.869986\n",
".. ... ... ... ... ...\n",
"995 -55.100100 34.979980 0.0 4.975266 22.643122\n",
"996 -55.075075 34.984985 0.0 4.975266 22.643122\n",
"997 -55.050050 34.989990 0.0 4.975266 22.643122\n",
"998 -55.025025 34.994995 0.0 4.975266 22.643122\n",
"999 -55.000000 35.000000 0.0 4.975266 22.643122\n",
"\n",
"[1000 rows x 5 columns]"
]
},
"execution_count": 15,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"particle_df_from_parcels = pd.concat(map(particle_to_pandas, pset))\n",
"\n",
"particle_df_from_parcels"
]
},
{
"cell_type": "code",
"execution_count": 16,
"id": "767b65a5-a0a2-4dac-8787-8a12bbc098fa",
"metadata": {},
"outputs": [],
"source": [
"import xoak"
]
},
{
"cell_type": "code",
"execution_count": 17,
"id": "b59d2d40-ec6c-4d02-bdb3-1eedc5ea82ab",
"metadata": {},
"outputs": [],
"source": [
"ds_T = xr.open_dataset(\"GYRE_5d/GYRE_5d_00010101_00011230_grid_T.nc\")\n",
"ds_T.xoak.set_index(['nav_lat', 'nav_lon'], 'sklearn_geo_balltree')"
]
},
{
"cell_type": "code",
"execution_count": 18,
"id": "b9d5a315-16b6-4049-8e13-3f52b447dad2",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>temperature</th>\n",
" <th>lon</th>\n",
" <th>lat</th>\n",
" <th>depth</th>\n",
" <th>time</th>\n",
" </tr>\n",
" <tr>\n",
" <th>index</th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>22.839523</td>\n",
" <td>-80.281670</td>\n",
" <td>30.348095</td>\n",
" <td>4.975266</td>\n",
" <td>0001-01-03 12:00:00</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>22.839523</td>\n",
" <td>-80.281670</td>\n",
" <td>30.348095</td>\n",
" <td>4.975266</td>\n",
" <td>0001-01-03 12:00:00</td>\n",
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <td>22.839523</td>\n",
" <td>-80.281670</td>\n",
" <td>30.348095</td>\n",
" <td>4.975266</td>\n",
" <td>0001-01-03 12:00:00</td>\n",
" </tr>\n",
" <tr>\n",
" <th>3</th>\n",
" <td>22.869986</td>\n",
" <td>-79.607620</td>\n",
" <td>29.674047</td>\n",
" <td>4.975266</td>\n",
" <td>0001-01-03 12:00:00</td>\n",
" </tr>\n",
" <tr>\n",
" <th>4</th>\n",
" <td>22.869986</td>\n",
" <td>-79.607620</td>\n",
" <td>29.674047</td>\n",
" <td>4.975266</td>\n",
" <td>0001-01-03 12:00:00</td>\n",
" </tr>\n",
" <tr>\n",
" <th>...</th>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" </tr>\n",
" <tr>\n",
" <th>995</th>\n",
" <td>22.643122</td>\n",
" <td>-55.341927</td>\n",
" <td>35.066425</td>\n",
" <td>4.975266</td>\n",
" <td>0001-01-03 12:00:00</td>\n",
" </tr>\n",
" <tr>\n",
" <th>996</th>\n",
" <td>22.643122</td>\n",
" <td>-55.341927</td>\n",
" <td>35.066425</td>\n",
" <td>4.975266</td>\n",
" <td>0001-01-03 12:00:00</td>\n",
" </tr>\n",
" <tr>\n",
" <th>997</th>\n",
" <td>22.643122</td>\n",
" <td>-55.341927</td>\n",
" <td>35.066425</td>\n",
" <td>4.975266</td>\n",
" <td>0001-01-03 12:00:00</td>\n",
" </tr>\n",
" <tr>\n",
" <th>998</th>\n",
" <td>22.643122</td>\n",
" <td>-55.341927</td>\n",
" <td>35.066425</td>\n",
" <td>4.975266</td>\n",
" <td>0001-01-03 12:00:00</td>\n",
" </tr>\n",
" <tr>\n",
" <th>999</th>\n",
" <td>22.643122</td>\n",
" <td>-55.341927</td>\n",
" <td>35.066425</td>\n",
" <td>4.975266</td>\n",
" <td>0001-01-03 12:00:00</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"<p>1000 rows × 5 columns</p>\n",
"</div>"
],
"text/plain": [
" temperature lon lat depth time\n",
"index \n",
"0 22.839523 -80.281670 30.348095 4.975266 0001-01-03 12:00:00\n",
"1 22.839523 -80.281670 30.348095 4.975266 0001-01-03 12:00:00\n",
"2 22.839523 -80.281670 30.348095 4.975266 0001-01-03 12:00:00\n",
"3 22.869986 -79.607620 29.674047 4.975266 0001-01-03 12:00:00\n",
"4 22.869986 -79.607620 29.674047 4.975266 0001-01-03 12:00:00\n",
"... ... ... ... ... ...\n",
"995 22.643122 -55.341927 35.066425 4.975266 0001-01-03 12:00:00\n",
"996 22.643122 -55.341927 35.066425 4.975266 0001-01-03 12:00:00\n",
"997 22.643122 -55.341927 35.066425 4.975266 0001-01-03 12:00:00\n",
"998 22.643122 -55.341927 35.066425 4.975266 0001-01-03 12:00:00\n",
"999 22.643122 -55.341927 35.066425 4.975266 0001-01-03 12:00:00\n",
"\n",
"[1000 rows x 5 columns]"
]
},
"execution_count": 18,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"particle_df_from_direct_selection = ds_T.xoak.sel(\n",
" nav_lon=particle_df_from_parcels.to_xarray()[\"lon\"],\n",
" nav_lat=particle_df_from_parcels.to_xarray()[\"lat\"],\n",
").sel(\n",
" deptht=particle_df_from_parcels.to_xarray().depth,\n",
" method=\"nearest\",\n",
").isel(\n",
" time_counter=0\n",
")[[\"votemper\"]]\n",
"\n",
"particle_df_from_direct_selection = (\n",
" particle_df_from_direct_selection\n",
" .to_dataframe()\n",
" .rename(columns={\n",
" \"votemper\": \"temperature\",\n",
" \"time_counter\": \"time\",\n",
" \"deptht\": \"depth\",\n",
" \"nav_lat\": \"lat\",\n",
" \"nav_lon\": \"lon\",\n",
" })\n",
")\n",
"\n",
"particle_df_from_direct_selection"
]
},
{
"cell_type": "code",
"execution_count": 19,
"id": "1d7f69fb-a93f-4977-9d50-ec492897ef55",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"count 1000.000000\n",
"mean -0.000057\n",
"std 0.004203\n",
"min -0.029823\n",
"25% 0.000000\n",
"50% 0.000000\n",
"75% 0.000000\n",
"max 0.029566\n",
"Name: temperature, dtype: float64"
]
},
"execution_count": 19,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"(\n",
" particle_df_from_direct_selection[\"temperature\"]\n",
" - particle_df_from_parcels[\"temperature\"]\n",
").describe()"
]
},
{
"cell_type": "code",
"execution_count": 20,
"id": "5b09112a-68bb-4d71-855d-f510cb239b95",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAjgAAAGwCAYAAACkfh/eAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABB00lEQVR4nO3de3gUVZo/8G/n1kkgtECGNNEAwRswUYFkhLCDMDsaEBjHUUcUDeKFFTMIIeuKgPuY8ULEQZafi8AjAs4sOPLsALvOPmyG6EoWNwlISLhGHDUShDQxXDqRS0KS+v3BpO1KX0jTXenUeb+f58nzNN2nm/NWddd565xTpyyapmkgIiIiUkhEuCtAREREFGpMcIiIiEg5THCIiIhIOUxwiIiISDlMcIiIiEg5THCIiIhIOUxwiIiISDlR4a5AOLS1teHEiRNISEiAxWIJd3WIiIioEzRNQ2NjI5KTkxER4b+PRmSCc+LECaSkpIS7GkRERHQVjh07huuuu85vGZEJTkJCAoDLG6hXr15hrg0RERF1RkNDA1JSUlztuD8iE5z2YalevXoxwSEiIjKZzkwv4SRjIiIiUg4THCIiIlIOExwiIiJSDhMcIiIiUg4THCIiIlIOExwiIiJSDhMcIiIiUg4THCIiIlJOlyQ4K1euRGpqKmJjY5Geno6dO3f6LV9cXIz09HTExsZi8ODBWL16te71LVu2ICMjA9dccw169OiB4cOH49/+7d+MDIGIiIhMxPAEZ9OmTcjNzcWiRYtQUVGBsWPH4u6770ZNTY3X8tXV1Zg0aRLGjh2LiooKLFy4EHPmzMHmzZtdZfr06YNFixahtLQU+/fvx+OPP47HH38cf/nLX4wOh4iIiEzAommaZuR/MGrUKIwcORKrVq1yPTd06FDce++9KCgo8Cg/f/58fPjhh6iqqnI9N2vWLOzbtw+lpaU+/5+RI0di8uTJeOWVVzxea2pqQlNTk+vf7feycDqdvFUDERGRSTQ0NMBms3Wq/Ta0B6e5uRnl5eXIysrSPZ+VlYWSkhKv7yktLfUoP2HCBOzZsweXLl3yKK9pGj7++GMcOXIEd9xxh9fPLCgogM1mc/3xTuJERERqMzTBqa+vR2trK5KSknTPJyUlweFweH2Pw+HwWr6lpQX19fWu55xOJ3r27ImYmBhMnjwZ//qv/4q77rrL62cuWLAATqfT9Xfs2LEgIyNSh6Zp+K6x6coFiYhMpEvuJt7xrp+apvm9E6i38h2fT0hIQGVlJb7//nt8/PHHyMvLw+DBgzF+/HiPz7NarbBarUFEQKSuf/rTfvyp/FvM+fmNyLvrpnBXh4goJAztwUlMTERkZKRHb01dXZ1HL007u93utXxUVBT69u3rei4iIgI33HADhg8fjn/8x3/EAw884HVOD109TdOw9tNqlB89E+6qkIH+VP4tAGDdp9VhrgkRUegYmuDExMQgPT0dRUVFuueLioowZswYr+/JzMz0KL99+3ZkZGQgOjra5/+laZpuIjEFb9sBB175r8O4f5X3+VIqutDcigVb9mPHkbpwV6XLtbS1hbsKREQhY/hl4nl5eXj33Xexbt06VFVVYd68eaipqcGsWbMAXJ4fM336dFf5WbNm4ejRo8jLy0NVVRXWrVuHtWvX4rnnnnOVKSgoQFFREb7++mt8/vnnWLZsGf7whz/g0UcfNTocUb767vtwV6HLrS7+Cn/cfQwz1n8W7qoQEVEQDJ+DM3XqVJw6dQovv/wyamtrkZaWhm3btmHgwIEAgNraWt2aOKmpqdi2bRvmzZuHt99+G8nJyXjrrbdw//33u8qcO3cOOTk5+PbbbxEXF4chQ4Zgw4YNmDp1qtHhkOKOn70Q7ioQEVEIdMkk45ycHOTk5Hh97b333vN4bty4cdi7d6/Pz3v11Vfx6quvhqp6REREpBjei4qIiIiUwwSHiIiIlMMEh4iIiJTDBIeIxGlt03CuqSXc1SAiAzHBISJxfvn2p/jxS39B/fdy1s568T8OIHvtLrS1GXp/ZQqj75ta8MfdNTjwrTPcVekWmOAQkTgHjzcAAP7nczkLOm4oq8HOv9ajvIYrk6vqneKvsGDLAfxixafhrkq3wASHiEiQViE9OK1tGkq+rMex0+fDXZUuUyMo1s5ggkNERMrZvPdbTHt3F36+rBjNLbwNiURMcIiISDnf/q03o7mlDRdbWsNcGwoHJjhERKQ0TcaonI4mMegOmOAQEZHa2NaLxASHiIiUpgnMcNiBwwSHSIcHBSJSAQ9lTHCIiEhxEk9cOAeHCQ6RjsUS7hoQUaixqZeJCQ4RESlNYm+GvIg9McEhIiJSjMCczgMTHPKJPxAiUoHEQ5nEK8c6YoJDRKQ46UM0AsMnMMEhPzjhlojMyj2pkdibwaSOCQ4RkfLY2JFETHCIiBQnMb/R9doI3ABMapngEOnwoEAqEjkHR3Z+Q2CCQ0RECpI+yVjivKOOmOAQueHEalIRmzp5JCZ1HTHBISJSnMTGTvxVVOGuQDfABIeIiJTjntRITPCICQ75wYMCkRok9mAIv4hK5MTyjpjgEBEpjm2dPNzlTHDID064JSKz0l9FJa+5FxiyByY45BN/IERkVu5JDY9lMjHBISJSnMQGXmLMOtLjBxMc8kPiEJX4gyIpSeQkY+G4z5ngEBEpT2LiLn4lY4Exd8QEh8iNxF4rIhVJX+iPmOAQ6fCsh1Qk8WstfaE/gSF7YIJDRKQ4iZdJS8d9zgSH/JD4++AQFalI4E+5wxAVScQEh4iIlCaxN0NexJ6Y4JBP7M0gUoPA9l2/0F8Y69GVpF851hETHPKJPxAiRfC3LAKvHNNjgkNEokgcrpCIvRnEBId8kjhExQMhqUji2bz+tywjfoEh+8UEh3xiY08qkvi9FhmzwHVwJM478ocJDpEbib1W0kg88EuMWSIOy+l1SYKzcuVKpKamIjY2Funp6di5c6ff8sXFxUhPT0dsbCwGDx6M1atX615fs2YNxo4di969e6N379648847sXv3biNDEEliY8+DApEauA4OGZ7gbNq0Cbm5uVi0aBEqKiowduxY3H333aipqfFavrq6GpMmTcLYsWNRUVGBhQsXYs6cOdi8ebOrzI4dO/Dwww/jk08+QWlpKQYMGICsrCwcP37c6HCIyOQkTjIWGbP7Yynh8yoqHcMTnGXLluHJJ5/EU089haFDh2L58uVISUnBqlWrvJZfvXo1BgwYgOXLl2Po0KF46qmn8MQTT2Dp0qWuMhs3bkROTg6GDx+OIUOGYM2aNWhra8PHH39sdDikOIm9VtJIPOyLjFlgYy9x3pE/hiY4zc3NKC8vR1ZWlu75rKwslJSUeH1PaWmpR/kJEyZgz549uHTpktf3nD9/HpcuXUKfPn28vt7U1ISGhgbdHxHJJPHAr2vsBcYvBYfl9AxNcOrr69Ha2oqkpCTd80lJSXA4HF7f43A4vJZvaWlBfX291/e88MILuPbaa3HnnXd6fb2goAA2m831l5KSchXREBGZn5TeDLA3Q7wumWRs6dDvr2max3NXKu/teQB444038Mc//hFbtmxBbGys189bsGABnE6n6+/YsWOBhiASDwqkIjkN/A9Exiyw10ofs5Cg/Ygy8sMTExMRGRnp0VtTV1fn0UvTzm63ey0fFRWFvn376p5funQpFi9ejI8++gi33nqrz3pYrVZYrdarjIKIVCLyuK/5eKwwzsEJY0W6CUN7cGJiYpCeno6ioiLd80VFRRgzZozX92RmZnqU3759OzIyMhAdHe167ne/+x1eeeUVFBYWIiMjI/SVJ5ETbnlQINXxK64uHr/0DB+iysvLw7vvvot169ahqqoK8+bNQ01NDWbNmgXg8vDR9OnTXeVnzZqFo0ePIi8vD1VVVVi3bh3Wrl2L5557zlXmjTfewIsvvoh169Zh0KBBcDgccDgc+P77740ORxT+WIjUIPGSafZmkKFDVAAwdepUnDp1Ci+//DJqa2uRlpaGbdu2YeDAgQCA2tpa3Zo4qamp2LZtG+bNm4e3334bycnJeOutt3D//fe7yqxcuRLNzc144IEHdP/XSy+9hPz8fKNDIoVJ7LWSRmJjx5hlkJjI+mN4ggMAOTk5yMnJ8frae++95/HcuHHjsHfvXp+f980334SoZuSPxMaeBwX1SZmP4U7XmyEkfomNvcR5R/7wXlTUKZyRT6qQ/lWWHr/aOCznjgkOkRuJvVakPomNHXsziAkOdYrEAySpSeJXWeBV4iInGXMlYz0mONQp/LGoSeLQo/SYxcQvsLHXzzuSErVvTHDIJ66KqT7uVnm4y0kKJjhEgkls7ETGLDBoib0Zup66MNaju2CCQz65T7iV8mMRchwUTfw+FhK/xMZe4qXx/jDBIZ8k3qxOGilntjoSQxZ4RZHExl4fp5Cg/WCCQ+RG2mXiPATKI6WxJ2KCQz7ph6hkHBV58FeflO+yO5ExC+zNkNhr5Q8THPKJQ1Tqk7hfpccsJX6Jjb3EeUf+MMEhciNviEreYVBexHpS4pfe2EtJ6vxhgkM+SWvsAR4USE0SL5kmYoJDncJjopok7leJDbzImN0fCwlf4tVy/jDBIZ/4YyEVSfwmS7wXle5WDUIyHIn33/KHCQ6RG4nDcqIJbASkNHy6xj6M9ehKEieT+8MEh3zSXSYu5MciJc520uIFZPZMSrxkmogJDvmkbwhIRVIaeHcSY5b4C5bYmyExefeHCQ6RGw5RCSC94ZMYs5DGnnNw9JjgkE/6ISoZvxYhYbpIi7cjieFLiVmX1AgJWvrvuSMmOOQTh6jUJ3G/irx8ONwVIAoDJjhEbjhEpT4pSY07DlHJIDF594cJDvnEq6jUJ2Xo0Z3+8mHZ8atMZGMvcN6RP0xwyCfdQYG/FSVJ361SGj6Jk085yTiMFekmmOBQp0g5QEgj8SAocuhCSqBEbpjgELmReOWYNBL3qsSkDgJ7M2TuZ9+Y4JBP0ufgiIhZQoz+iNjJelISd4mNPe8ar8cEh3ySeICQRuLQo/uBX0r0Ivez+2M29iIxwSFyo+u1Cl81yEAS2zqJMUskMXn3hwkO+SRxPop+iEr9mAWE6JfE+KXELLGxF3lpvB9McMgnDlGpT/p+lZDEdiRluEriDdR513g9JjhEbjhEpT6Jibu4yfOQuQ4O6THBIZ8kLholrSFgD4YMEmOWiENUekxwyCfpZ0ASYlY/Qk/SklhAaMzuj4XEDIHzjvxhgkM+cThXFjGNgBuBIYuJWTfJWEjQIpM6P5jgELmRtrihyPko7o8l7GTIjNmdvIgJYIJD/gjs7pR27JcwDNeRyAZeYMwSSVvm4kqY4JBP7O4URsg+FhKmT1Lil9jY6y4MCWM9ugsmOERupA1RQTdEJSFgPRH7GDLn00ls7CVOJveHCQ75JPEqKmkxqx+hJ2n7GJDZ2LGxJyY45JPEdXCkkdkIiAnUjXtvhsT4ZZCYvPvDBIfIjbghKjfCwgUgbx8DcmKWeNsCgSH7xQSHfBJ5CbGwmPW9dBIilrePAaExC+yBlniDUX+Y4JBP4tfOEBCzgBA9CAxZZswCkzrSY4JDRADkNAIS5x1JjFk67mcmOOSHxIOitLuJS1/rSOJETCkxS/xuc5KxXpckOCtXrkRqaipiY2ORnp6OnTt3+i1fXFyM9PR0xMbGYvDgwVi9erXu9UOHDuH+++/HoEGDYLFYsHz5cgNrL5fEH4i0pE7CMFxH4udmCIlZ4hpPEr/b/hie4GzatAm5ublYtGgRKioqMHbsWNx9992oqanxWr66uhqTJk3C2LFjUVFRgYULF2LOnDnYvHmzq8z58+cxePBgvP7667Db7UaHIJewxp5kkPhdFhgyG3syPsFZtmwZnnzySTz11FMYOnQoli9fjpSUFKxatcpr+dWrV2PAgAFYvnw5hg4diqeeegpPPPEEli5d6irzk5/8BL/73e/w0EMPwWq1Gh0CCeI+RCWhVeh44JfQoyMgRL+Eh680TqzWMzTBaW5uRnl5ObKysnTPZ2VloaSkxOt7SktLPcpPmDABe/bswaVLl66qHk1NTWhoaND90ZXpxrCF/Fw4hi2LhIQO6JDUCYxZRsS88rUjQxOc+vp6tLa2IikpSfd8UlISHA6H1/c4HA6v5VtaWlBfX39V9SgoKIDNZnP9paSkXNXnSCNy3F44CftZ4tCFyPsyuT+WsqNJp0smGVt0/f6Xv2wdn7tSeW/Pd9aCBQvgdDpdf8eOHbuqzyH1SVvJ2GOIKjzV6FIS9qsHgTFLTGq40J9elJEfnpiYiMjISI/emrq6Oo9emnZ2u91r+aioKPTt2/eq6mG1WjlX5yqI7OIVFrP0YTiJ0Qts98XgrRr0DO3BiYmJQXp6OoqKinTPFxUVYcyYMV7fk5mZ6VF++/btyMjIQHR0tGF1JU/Su3gZs/qkhCvxtyxw2pHIS+P9MXyIKi8vD++++y7WrVuHqqoqzJs3DzU1NZg1axaAy8NH06dPd5WfNWsWjh49iry8PFRVVWHdunVYu3YtnnvuOVeZ5uZmVFZWorKyEs3NzTh+/DgqKyvx5ZdfGh0OkVKkD1FJaQSk9UwCMvcz6Rk6RAUAU6dOxalTp/Dyyy+jtrYWaWlp2LZtGwYOHAgAqK2t1a2Jk5qaim3btmHevHl4++23kZycjLfeegv333+/q8yJEycwYsQI17+XLl2KpUuXYty4cdixY4fRIYkh8aDoTkLMEmLsSGJjJzNmt8dCwpcYsz+GJzgAkJOTg5ycHK+vvffeex7PjRs3Dnv37vX5eYMGDRLTzRpO4q82ERKzO2kxS4sXkBmzFLzyVY/3oiKf9D8Q/lpU1PFEQcKZvsSeSYkxQ2Bjr/l4LBUTHCJ3wsbtO0YooSHQ5+0CAobMxo6NPTHBoU4R0g6Iu8xSyn51J3F4Wz90ISN+/Y1zGbNETHDIJ4mLRkmMWTLuY1KJxBWr/WGCQyRahzk4Ao6KEq80kRmzvMZe34MTvnp0F0xwyCeZB0W3xwJilhBjRxJjFtPCu9HfYDRs1aAwYoJDPklcKEtizO6kxSwtXkBOzBJ/y7zyVY8JDpFgEq+ico9aRrxc30ki7mcmOOSHxIMih6jUJ3FNGJExuz8WEjQvktBjgkM+SZywJv0AISFmCTF2JOX36076b5mY4JAf+jl6PESoqON+lbZ2hrBwATBmlUnstfKHCQ6RG/0BQv0jhIAQPYicfKp7LCNmiSR+t/1hgkM+SRyigsSY3UgIWRM4IUXiTRglNvYS5036wwSH/JD3C5FyIGzX8SAo4aAoIEQPMmNmYy8dExwiwTwSOmENgbBwRRHYUScyZn+Y4JBPEoeoJMYsjcQbEkqMWSJpcwivhAkO+SRyDFtYzB5DVBJiFjl0IS9m/aK+QoImHSY4ROQioh2QECOJXAdHxO83AExwyCeJZ7rSYpYQoz9Swpc4N0PmmjCyjl9XwgSHfBJ5UBQYszsJMUts+CTGLJG0IfYrYYJDPnHCmvokrmQssREQsFs9cWK1eExwiNxIS+oEhEgdiEnqfDxWGXvq9JjgkE8Sh2skxuxOQszS5lkBQmMWuXqzvJj9YYJDPkk8KEqbpNcxRBExC4ixI+mJuxQSe638YYJDJJiEYbiO5EXMmCXGT0xwyB/N5z+UpQmM2Z2U+RntJCZ4UrqwJK7eLDFmf5jgkE8SJ6xJi9kjRAkxi1wATmDMYiL9gcT97A8THCLBPG/VoD5pSWxHEmOWQngHtAcmOOSTxLMBiTGLI3DHSkxqeONcYoJDPkk805UXc8eF/sJUjTAROYwhJGaJCzpCYsx+MMEhnzhhTX28m3gYK9KFJMbsTkrM8k7Q/GOCQ+RG5FmfMBLXhJEYs0QcYtdjgkM+SVxHQtoZkMSF/kgGNvbEBId8krjst7SYRV5FJXDyqciY3R8zZpGY4BCRKJqff6lK3xsrJGaBw80SY/aHCQ75JPGg6E5CzB0nj0uYTC6tlw7osF+FxCyR9MnkHTHBId8ETsKR1pXPOTikKjb2xASHSDCJB36J8xQEnquI2bfu2FGnxwSHfNKdAYWxHl1JwrCUdCLnKQhc00qfyMqLWWSG1wETHPJJ2nANIC/mjg28hJghcOhCTCJH5IYJDvkk8kxXGoErGUsnJqkTdrICoMOtGogJDpEbJnXqkzhPQWLMED7ELiap84MJDvkk8cciLWaJV1GJn2QsJWaBPTi8f6AeExzySeJZn7SYJa5kTEQyMMEhIlEkDkOKjFn3WGLMxASHfJJ+maWEmD2vohIWs/rhApA39ApwxWopMfvTJQnOypUrkZqaitjYWKSnp2Pnzp1+yxcXFyM9PR2xsbEYPHgwVq9e7VFm8+bNGDZsGKxWK4YNG4atW7caVX2xpA3XABB3FYLEISqJ32uJjZ3E3gyJMftjeIKzadMm5ObmYtGiRaioqMDYsWNx9913o6amxmv56upqTJo0CWPHjkVFRQUWLlyIOXPmYPPmza4ypaWlmDp1KrKzs7Fv3z5kZ2fjwQcfxK5du4wOh4iIiEzAohncJz1q1CiMHDkSq1atcj03dOhQ3HvvvSgoKPAoP3/+fHz44YeoqqpyPTdr1izs27cPpaWlAICpU6eioaEB//3f/+0qM3HiRPTu3Rt//OMfPT6zqakJTU1Nrn83NDQgJSUFTqcTvXr1CkmcAHDs9Hk8tm53p8t/XX8OiT1j0Cs2OmR1CCVHw0Wcb24FACT1sqJHTFTQn/l1/Tkk22IRGx0Z9GcZ4fjZC2hqaQOAkNQzMsKCiy2tiI7onqPB55tb4Wi46Pr3tdfEwRoVXF2/rj8HABic2COozzFKY1MLvmu8fDzoaY1CvwRr0J954VIrromPQdOl1qA/ywgNFy+h/vtmAEBCbBR+1DP4mJta2pAQG4Xmv/1eupvqU+dcPVfXxEejT3xMUJ/Xqmk4eup8t/1eAz/89gCgT48YXBMXXNsSExWBxostV31MuK5PPP7wxO1B1aGjhoYG2Gy2TrXfwbdYfjQ3N6O8vBwvvPCC7vmsrCyUlJR4fU9paSmysrJ0z02YMAFr167FpUuXEB0djdLSUsybN8+jzPLly71+ZkFBAX77299efSCddKm1TfcF64z675tdB57u7GRDE4CmK5brjBPOi1cu1A2YpZ6hdPzshZB9VqC/hXD4vqkF3ze1hOSzak3yfWm82ILGi6GJ2SzOnr+Es+cvheSzzPC9BoDT55px+lyY2xZLeP97QxOc+vp6tLa2IikpSfd8UlISHA6H1/c4HA6v5VtaWlBfX4/+/fv7LOPrMxcsWIC8vDzXv9t7cEIt+Zo4/PuszE6VfXHrQRw52QgAWPyrW3BjUs+Q1ycUbHHRaLzYgrYQdPTlflDpakD/30PDkXxNXNCfaYRr4qLhvHAp6DHsnX+tx1sf/9X177WPZaBXkGdURoiwADf8KAFf1DUG/VmFBx1Y+2k1AGDSLXY8/nepQX+mESIjLIiPiQxJQz9j3W6c+1tP56/Tr8ODPwn9sSUUYiIjEB0ZgXPNwcf869WlrsczxgzC5Fv7B/2ZRhjYJx51jU24EGTP2vcXW/D4e58BAOJjIvH7EPdKhFJqYg+ccOuJvlpr/vdrbD98EsDl3uz/9/CIgD8jNiq8PfWGJjjtLBZ9GqdpmsdzVyrf8flAPtNqtcJqDb5L9kpioyPxk0F9OlU2IfaHTT8suReGp1xjUK26j/iYH77st153DVK7cVdvKJzo0BsyYkBv9OkRXDe5kTr73fXn89oG1+OkXrEh+czuLjLih+POtb3jRMTsbmDf+G4dc79esUF/htOt9yc6MqJbxwsAiSEYgvyw8oTrcWxM59u27sTQiQGJiYmIjIz06Fmpq6vz6IFpZ7fbvZaPiopC3759/Zbx9ZndkXsuFuZevC4jLeaOCbeEmN13skVGxPoTLyExuxMRsYgg9VQ4Xhua4MTExCA9PR1FRUW654uKijBmzBiv78nMzPQov337dmRkZCA6OtpvGV+fSURkFF1DYNaWIED6mNUPWkCIHtxDNus+NnyIKi8vD9nZ2cjIyEBmZibeeecd1NTUYNasWQAuz485fvw4/vCHPwC4fMXUihUrkJeXh5kzZ6K0tBRr167VXR01d+5c3HHHHViyZAl++ctf4j//8z/x0Ucf4dNPPzU6nJBxP9Mz6XcnYNJi7hiitJglxEtE3ZfhCc7UqVNx6tQpvPzyy6itrUVaWhq2bduGgQMHAgBqa2t1a+KkpqZi27ZtmDdvHt5++20kJyfjrbfewv333+8qM2bMGHzwwQd48cUX8c///M+4/vrrsWnTJowaNcrocEJH1/0noyWwCIu5YwMvLWb1o73M4uOxyiz4YSE5CYmsgBA96IdezalLJhnn5OQgJyfH62vvvfeex3Pjxo3D3r17/X7mAw88gAceeCAU1aMwkHFQ9MhwSEH6ix/CWJEuZLFYXMsjSwjZrEM0oWLW8Lvn6mMCsCufVCRtGJJkkPhVVqHHnQlOmEg8+Es7C/IYohIQvrTJp4AakzEDpYtSSMxkPkxwwkTima60XiuPScZhqUXXEjkfRUqgbqTNtRK5jxVoo5jgEBEFxfwNQaBUaPwCYdYhmmCosF+Z4ISJCuObgZI2fOE5RCUsZvXDJSEE/HQ9qDD0ygQnTMQvDha+anQheSsZ687sRUQs82RF4jIXZD5McIiIgiBtbhkgM2ZpVDghZYITJtLGsAF5MUu8igrsmRRBWsxSvsvuVFjfiQlOmEjs1pYWs+dVVLJiVj9akkLCb7cjFXrpmOBQl1HhBxMIj7uJC4hZImk9k4DMmMl8mOB0AzxAkCpU6NYOlLSeSUBezFK+yzoK7GMmOGGiwo3MAiYsZgkxdqQfopKxBaT1TAIdVzIOVy26joAQPajQS8cEh7qMuIOixEnGRAoy6zowwVBhIjkTnDCReNZH6pO5vpOQQN2I7IEm02GCEyYWad0Z4Li9tJjVj9aTlGRHhVVuA6F+hJ70TZQ5twATHOoy0nqtOiY04mKWEDCEJnXCYpbZM+n2OHzVCAoTnDCR1tgTEZmVey+VpoWxIl2Ik4zpqkkcwxYXs8cQlfpUOOsLlMize/fHQmIm82GCQ2EhcdxeQszupISrv8GoDBLXO5JGhZMVJjhhIm2SXkcSIvZYyThM9SCi0JFyuFahjWKCEyYqZMdEHemHIWV8s/VDVAJjFrKf20mZg+POrHuYCU7YyO7ilRCz5xBVWKrRpSTOzWDMYasGGUmBHcsEh8JCwlmf50rG6scskbjJ8+B3WQIVklgmOGEiuYuX1MWhV1KdWRv7QKnQRjHBCRMVsuOgCIjZrAeFYKiwdkagVFjxNVAqTEC9WhLn4Jj1UMYEh8JCwjFR4vooIgnstWJPnfpUOEFjghMmbPDUZ/HxWGUiryhyfywjZEi/SEICFZJYJjhhIrEr353AkEXgfiXVSTleq5C4M8GhsBBxdi+xN0PgsJz4tX+ExNxO4hwcs+5jJjhhIrErXxqJS/hLpMKZbqAkxiyNCvuVCU6YqDC+GQyJMcsgvTdDBokxS6PC/caY4FBYmPUHEwiZwzXeHxOpQuL32qwxM8EJE04yVj9o/VVU6scrlcTfssSY23EOjnkwwQkXwZP0SF3iL40XErU+qZERszQqJK5McCgsVPjxXIlF4EQFFcbtgyIkZk4yVp8KvXRMcMKEBwj1CcxviJQn5XitQpxMcMJE4h2IpZGYxEqcdyTxtywx5nYi5+CY9ADGBCdMpA9hm/T3Qlcg8SoqyTeeBGTGLIEKe5UJDoWFhLN7iZNPiUgNKgyxM8EJEzZ+Eph/kl6gpMTpToWGIFASY24n5TvOScZ01STOz3AnMWYJ9AdFGTtZ5LCcwJjbiZyDE+4KXCUmOBQWZv3BBELkWa7AmFU40w2UxJilUWG/MsEJE8lXIQAyzu6lTz4lIjWY9fjFBKcbMOuXh6gjiUOvEufTSYy5nZzvtflPwpngUFiY9QcTCBUOEIESGbP7YyFBS17mQuQcHJPuY0MTnDNnziA7Oxs2mw02mw3Z2dk4e/as3/domob8/HwkJycjLi4O48ePx6FDh3Rl3nnnHYwfPx69evWCxWK54md2dyb97gTFrD+YQEhuBESR8GXuQGIiK43Fz7/MwtAEZ9q0aaisrERhYSEKCwtRWVmJ7Oxsv+954403sGzZMqxYsQKfffYZ7HY77rrrLjQ2NrrKnD9/HhMnTsTChQuNrH6XEXh8FEHiJGPOOyJSgwo/3yijPriqqgqFhYUoKyvDqFGjAABr1qxBZmYmjhw5gptvvtnjPZqmYfny5Vi0aBHuu+8+AMDvf/97JCUl4f3338fTTz8NAMjNzQUA7Nixw6jqk8HY+KlJ4uXDEpM6iTG3kxKuCkOvhvXglJaWwmazuZIbABg9ejRsNhtKSkq8vqe6uhoOhwNZWVmu56xWK8aNG+fzPZ3R1NSEhoYG3V93Im2SnhQS14SRSGJPncTlANqJnIMT7gpcJcMSHIfDgX79+nk8369fPzgcDp/vAYCkpCTd80lJST7f0xkFBQWueUA2mw0pKSlX/VmGMOu3h/yS2Zshb26GCme6gZIYszS6eVYm3ccBJzj5+fmwWCx+//bs2QPA+1mrpmlXPJvt+Hpn3uPPggUL4HQ6XX/Hjh276s8iIt8sbPmIlKDCzzfgOTizZ8/GQw895LfMoEGDsH//fpw8edLjte+++86jh6ad3W4HcLknp3///q7n6+rqfL6nM6xWK6xW61W/32gqfJHIP+5idemvKJKxpyXG3E7K8Vp3rmLSfRxwgpOYmIjExMQrlsvMzITT6cTu3btx++23AwB27doFp9OJMWPGeH1Pamoq7HY7ioqKMGLECABAc3MziouLsWTJkkCrahrm/OrQleiHqGTsZYuPxyqT2GklMeZ2IufgmHQfGzYHZ+jQoZg4cSJmzpyJsrIylJWVYebMmZgyZYruCqohQ4Zg69atAC43Arm5uVi8eDG2bt2KgwcPYsaMGYiPj8e0adNc73E4HKisrMSXX34JADhw4AAqKytx+vRpo8Ihos4QOO9IIpETq6VRYA6OYZeJA8DGjRsxZ84c11VR99xzD1asWKErc+TIETidTte/n3/+eVy4cAE5OTk4c+YMRo0ahe3btyMhIcFVZvXq1fjtb3/r+vcdd9wBAFi/fj1mzJhhYETGkHJ2L43MCbfyhi7Y2JOKVPguG5rg9OnTBxs2bPBbRuvQ32exWJCfn4/8/Hyf77nS62ajwheJPEm8ikoiiXfW1iWvQmJuJ2YfK3C/Md6LqhuQ8oORRr9fZexkkUmd9P0sJOZ2EufgmHUXM8EhopCROMmYSEUqDLEzwekGpJ0BSSFx6EIiyVcUAfJilhKvCnEywekGVPgikSeJk09VWP00UOL3cxjrQcZR4X5jTHCIKGQkzs2QeM8xFRq/qyVxDo5Z9zATHCKDSB+6ICLzUuGCASY43YBZvzzkn8zeDF//UJfMISrvjyWQEq8KxywmOEQUMmzsw1ePriRxP4ujwD5mgtMNqJApkzfyJtxKJPFqOYkxtxM5B8ekO5kJTjdg0u8OXYHMs1x5E26JVKTCmlZMcIgMIvFKE4lJnci5VhInW/2NkJ+yEscsJjjdgPm/RkSXib9yTEjM4vezACrksExwugEVMmXyxP0qg8hF7yTG/Dci5+CYdC8zwSEyiMSzXIkrGROpSIWrA5ngdAMm/e4QedBPTJTxzRY518r9sZCY20kJV4U4meB0Ayp8kciTCmdAdGWcWE0q4t3EKSSknQFJoT9AyNjHTOpkkDj82k7kHByT7mMmOEQGkdjYS0nk3Elu7EldKix/wASHiCgI+quozNkQBEpizO2YxJoHExyiLiDlmKjvtZIRtcQeHIkxS6PCFZFMcIgMIrGxdycvYpKAc3DMgwkOEVEQJF5RZNYGjzrP4udfZsEEh8ggEle4lTixWrd3hcQs+W7iUuJVIU4mOEQGUeFeLoESf2m8kJghMWZhVEhimeAQUcjI7MEhSUTOwQl3Ba4SExwig0icmyGRxCuKJMYsjQonK0xwiAyi7+I16REiQBKTOvExSwn6b6TEq0KYTHCIDCKy4VNg3D5QIhNZgXOtpFFhbhkTHCIiok4SOQfHnPkNExwio0icp6CPU0bQEodrJMYsj/mXuWCCQ2QUBbp4AyU9qRMSssiY20n8XpsVExwioiBw3lEYK0KG0Z+smHMnM8EhMojERkDymT3JIHEOjlkxwSEyiJSkRk/eFUUSl6yWONdKGt5NnIhIOInzjtxJi1lKvCqEyQSHyCAqjGEHikNUMki8kaw0XAeHiHySktS4k9ibwcZeFolzcMz6W2aCQ2QQkx4TgqLCuH2gRPbUuT8WErM0KvTGMsEhIgqCCg1BoCTG3E5KPmfWYSl3THCIDCLlQOhOd2avwAGSSCwFVqtmgkNkEIkNvEWf4Yggct6R+2MhMbeTOQfHnDuZCQ4RURD0k4zN2RAESmLM0qhwrsIEh8goZj0qBEG3enMY69GVJPZmSIy5nZR4zdpr444JDpFBFDg+BEx/l2mBG4BIESosVs0Eh4goGCY9+AdDcu4qcg6OSb/kTHCIDGLOQ0LoSIlf4k1VITJmWfS9seGrRzAMTXDOnDmD7Oxs2Gw22Gw2ZGdn4+zZs37fo2ka8vPzkZycjLi4OIwfPx6HDh1yvX769Gk8++yzuPnmmxEfH48BAwZgzpw5cDqdRoZCFDCJQzQqHBQDpcKS9oGSPBQpJVwVvsuGJjjTpk1DZWUlCgsLUVhYiMrKSmRnZ/t9zxtvvIFly5ZhxYoV+Oyzz2C323HXXXehsbERAHDixAmcOHECS5cuxYEDB/Dee++hsLAQTz75pJGhEAXM/IcHIpJKhcUco4z64KqqKhQWFqKsrAyjRo0CAKxZswaZmZk4cuQIbr75Zo/3aJqG5cuXY9GiRbjvvvsAAL///e+RlJSE999/H08//TTS0tKwefNm13uuv/56vPbaa3j00UfR0tKCqCjPkJqamtDU1OT6d0NDQ6jDJSLIvHxY4hVFKlxCTJ1n1u+1YT04paWlsNlsruQGAEaPHg2bzYaSkhKv76murobD4UBWVpbrOavVinHjxvl8DwA4nU706tXLa3IDAAUFBa5hMpvNhpSUlKuMiqjzzHpQCIbIxl76sJyQmKVRYVVywxIch8OBfv36eTzfr18/OBwOn+8BgKSkJN3zSUlJPt9z6tQpvPLKK3j66ad91mXBggVwOp2uv2PHjnU2DKKrZtaDQjBU6NamK9Ovd8Q9rSQFdmvACU5+fj4sFovfvz179gDwPvlM07QrTkrr+Lqv9zQ0NGDy5MkYNmwYXnrpJZ+fZ7Va0atXL90fEVEosLEn1Zm1ly7gOTizZ8/GQw895LfMoEGDsH//fpw8edLjte+++86jh6ad3W4HcLknp3///q7n6+rqPN7T2NiIiRMnomfPnti6dSuio6MDDYXIUGY9KARD18ALiV/icI3EmKVRYVXygBOcxMREJCYmXrFcZmYmnE4ndu/ejdtvvx0AsGvXLjidTowZM8bre1JTU2G321FUVIQRI0YAAJqbm1FcXIwlS5a4yjU0NGDChAmwWq348MMPERsbG2gYRGQA8ZdMh68aXUpizNLob5xrzr1s2BycoUOHYuLEiZg5cybKyspQVlaGmTNnYsqUKborqIYMGYKtW7cCuDw0lZubi8WLF2Pr1q04ePAgZsyYgfj4eEybNg3A5Z6brKwsnDt3DmvXrkVDQwMcDgccDgdaW1uNCoeIiP5GYk+dNCrsVsMuEweAjRs3Ys6cOa6rou655x6sWLFCV+bIkSO6Rfqef/55XLhwATk5OThz5gxGjRqF7du3IyEhAQBQXl6OXbt2AQBuuOEG3WdVV1dj0KBBBkZE1HkmPekJigInfVeBq/qS2sz6tTY0wenTpw82bNjgt4zW4cYeFosF+fn5yM/P91p+/PjxHu8h6o6kDNHoCDyx1yc1QqIWOBQpjW5NK5PuYt6LisggZj0oUGAk9lpJjFkaFebTMcEhopDR33jSnAdFIlIjiWWCQ2QQkx4TgiLx8mGJVxTpb8lB1D0xwSEyCHswZJDeayUxZglUSNyZ4BBRyEi8CaMKDUGg3C/0kBKzPJxkTEQ+mPSYEBSezROpQT/cbM7fNRMcIoOY9JgQFIEhKzEZM1AqXEJM6mOCQ2QQs571UGD0E27l7XOJMUugwl5lgkNEISM9p5MSv26xVSExS6NCLx0THCIKGZ7NE6lBf8GAOX/XTHCIiIJg1rPbYKhwdk/qY4JDRKEjsLGzKHA5bTAEhiyCCot2MsEhopAx64GQrh4n06tJl7iHsR7BYIJDRBQEFdYLIeqIPThERG5MehwMisTVm91JjJnMgQkOEYWMxB4MFc50gyExZml4FRUREYlj1saP/FMhcWeCQ0QhY9LjYFCkr2RMalLhu8wEh4hCxqxnesGQeC8qdxJjJnNggkNEIaPCWV/ALF4fKk13qwZSkgpXBzLBISIiIh2LAok7ExwiChmTnugFxSKwC4e3alCfCit0M8EhIgqC/kzXpC1BECTGTObABIeIKAgSJxm7z8GRErM0HKIiInLDxo5IDfrE3Zw/bCY4RBQyEocrVDjTDZR+7R9SERf6IyIiF7Oe6QZDYsxkDkxwiChkpLd1EsOXGLMM5u+lY4JDRCFj1gMhEelZFJg9zwSHiChETNoOBEVizBJYfDw2EyY4RBQy0udjyJxkLS9mMgcmOEQUMmzqiNSgwmrVTHCIKGTMeiAMGenxkzIi3X7MESb9YUeFuwJERKowaTtAARhiTwh3FbrEdb3jcHeaHbXOixh7Y2K4q3NVmOAQUchIn48hO3q1/ffcsfi3sqOY+/Mbw12VLhERYcGqR9PDXY2gMMEhIiK6gqH9e2Hxr24JdzUoAJyDQ0QUhGviYlyPoyNlHFLbh2l6x0eHuSZEvrEHh4goCFNvT0FLm4ab7T0RGx0Z7up0iX+4YzB+nNwLNyXJmI9C5sQEh4goCL1io/HM+OvDXY0uFR0ZgfE39wt3NYj8ktGfSkRERKIwwSEiIiLlMMEhIiIi5TDBISIiIuUwwSHqAkKuHiYi6jZ42KUuc02cvDUzfnFbMiIswLTbB4a7Kl1OC3cFiEg0QxOcM2fOIDs7GzabDTabDdnZ2Th79qzf92iahvz8fCQnJyMuLg7jx4/HoUOHdGWefvppXH/99YiLi8OPfvQj/PKXv8Tnn39uYCQUCi/fm4aMgb2x8pGR4a5Kl3nroeE4+NsJmDZqQLirQkQkiqEJzrRp01BZWYnCwkIUFhaisrIS2dnZft/zxhtvYNmyZVixYgU+++wz2O123HXXXWhsbHSVSU9Px/r161FVVYW//OUv0DQNWVlZaG1tNTIcCtK118ThT8+MwaRb+oe7Kl3GYrEgPkbmclP9bbHhrgIRCWbRNM2QnuSqqioMGzYMZWVlGDVqFACgrKwMmZmZ+Pzzz3HzzTd7vEfTNCQnJyM3Nxfz588HADQ1NSEpKQlLlizB008/7fX/2r9/P2677TZ8+eWXuP76Ky+41dDQAJvNBqfTiV69egUR5dWb/6f92LTnGADgm9cnh6UOREbYW3MGzvOX8LMhXAiOiEIrkPbbsB6c0tJS2Gw2V3IDAKNHj4bNZkNJSYnX91RXV8PhcCArK8v1nNVqxbhx43y+59y5c1i/fj1SU1ORkpLitUxTUxMaGhp0f+GWl3UTrv9RD7w4eWi4q0IUUiMH9GZyQ0RhZ1iC43A40K+f50GuX79+cDgcPt8DAElJSbrnk5KSPN6zcuVK9OzZEz179kRhYSGKiooQExMDbwoKClzzgGw2m89EqCsl9YrFx/84Hk+NHRzuqhARESkn4AQnPz8fFovF79+ePXsAXJ5/0JGmaV6fd9fxdW/veeSRR1BRUYHi4mLceOONePDBB3Hx4kWvn7dgwQI4nU7X37FjxwIJmYiIiEwm4NmPs2fPxkMPPeS3zKBBg7B//36cPHnS47XvvvvOo4emnd1uB3C5J6d//x8motbV1Xm8p7035sYbb8To0aPRu3dvbN26FQ8//LDH51qtVlit1ivGRkRERGoIOMFJTExEYmLiFctlZmbC6XRi9+7duP322wEAu3btgtPpxJgxY7y+JzU1FXa7HUVFRRgxYgQAoLm5GcXFxViyZInf/0/TNDQ1NQUYDREREanIsDk4Q4cOxcSJEzFz5kyUlZWhrKwMM2fOxJQpU3RXUA0ZMgRbt24FcHloKjc3F4sXL8bWrVtx8OBBzJgxA/Hx8Zg2bRoA4Ouvv0ZBQQHKy8tRU1OD0tJSPPjgg4iLi8OkSZOMCoeIiIhMxNAFOjZu3Ig5c+a4roq65557sGLFCl2ZI0eOwOl0uv79/PPP48KFC8jJycGZM2cwatQobN++HQkJCQCA2NhY7Ny5E8uXL8eZM2eQlJSEO+64AyUlJV4nNRMREZE8hq2D0511h3VwiIiIKDDdYh0cIiIionBhgkNERETKYYJDREREymGCQ0RERMphgkNERETKYYJDREREymGCQ0RERMphgkNERETKMXQl4+6qfW3DhoaGMNeEiIiIOqu93e7MGsUiE5zGxkYAQEpKSphrQkRERIFqbGyEzWbzW0bkrRra2tpw4sQJJCQkwGKxhPSzGxoakJKSgmPHjvE2EAbidu4a3M5dh9u6a3A7dw2jtrOmaWhsbERycjIiIvzPshHZgxMREYHrrrvO0P+jV69e/PF0AW7nrsHt3HW4rbsGt3PXMGI7X6nnph0nGRMREZFymOAQERGRcpjghJjVasVLL70Eq9Ua7qoojdu5a3A7dx1u667B7dw1usN2FjnJmIiIiNTGHhwiIiJSDhMcIiIiUg4THCIiIlIOExwiIiJSDhOcEFq5ciVSU1MRGxuL9PR07Ny5M9xVMpWCggL85Cc/QUJCAvr164d7770XR44c0ZXRNA35+flITk5GXFwcxo8fj0OHDunKNDU14dlnn0ViYiJ69OiBe+65B99++21XhmIqBQUFsFgsyM3NdT3H7Rwax48fx6OPPoq+ffsiPj4ew4cPR3l5uet1bufQaGlpwYsvvojU1FTExcVh8ODBePnll9HW1uYqw20duP/93//FL37xCyQnJ8NiseA//uM/dK+HapueOXMG2dnZsNlssNlsyM7OxtmzZ4MPQKOQ+OCDD7To6GhtzZo12uHDh7W5c+dqPXr00I4ePRruqpnGhAkTtPXr12sHDx7UKisrtcmTJ2sDBgzQvv/+e1eZ119/XUtISNA2b96sHThwQJs6darWv39/raGhwVVm1qxZ2rXXXqsVFRVpe/fu1X72s59pt912m9bS0hKOsLq13bt3a4MGDdJuvfVWbe7cua7nuZ2Dd/r0aW3gwIHajBkztF27dmnV1dXaRx99pH355ZeuMtzOofHqq69qffv21f7rv/5Lq66u1v793/9d69mzp7Z8+XJXGW7rwG3btk1btGiRtnnzZg2AtnXrVt3rodqmEydO1NLS0rSSkhKtpKRES0tL06ZMmRJ0/ZnghMjtt9+uzZo1S/fckCFDtBdeeCFMNTK/uro6DYBWXFysaZqmtbW1aXa7XXv99dddZS5evKjZbDZt9erVmqZp2tmzZ7Xo6Gjtgw8+cJU5fvy4FhERoRUWFnZtAN1cY2OjduONN2pFRUXauHHjXAkOt3NozJ8/X/vpT3/q83Vu59CZPHmy9sQTT+ieu++++7RHH31U0zRu61DomOCEapsePnxYA6CVlZW5ypSWlmoAtM8//zyoOnOIKgSam5tRXl6OrKws3fNZWVkoKSkJU63Mz+l0AgD69OkDAKiurobD4dBtZ6vVinHjxrm2c3l5OS5duqQrk5ycjLS0NO6LDn7zm99g8uTJuPPOO3XPczuHxocffoiMjAz8+te/Rr9+/TBixAisWbPG9Tq3c+j89Kc/xccff4wvvvgCALBv3z58+umnmDRpEgBuayOEapuWlpbCZrNh1KhRrjKjR4+GzWYLeruLvNlmqNXX16O1tRVJSUm655OSkuBwOMJUK3PTNA15eXn46U9/irS0NABwbUtv2/no0aOuMjExMejdu7dHGe6LH3zwwQfYu3cvPvvsM4/XuJ1D4+uvv8aqVauQl5eHhQsXYvfu3ZgzZw6sViumT5/O7RxC8+fPh9PpxJAhQxAZGYnW1la89tprePjhhwHwO22EUG1Th8OBfv36eXx+v379gt7uTHBCyGKx6P6taZrHc9Q5s2fPxv79+/Hpp596vHY125n74gfHjh3D3LlzsX37dsTGxvosx+0cnLa2NmRkZGDx4sUAgBEjRuDQoUNYtWoVpk+f7irH7Ry8TZs2YcOGDXj//ffx4x//GJWVlcjNzUVycjIee+wxVzlu69ALxTb1Vj4U251DVCGQmJiIyMhIj2yzrq7OI7ulK3v22Wfx4Ycf4pNPPsF1113net5utwOA3+1st9vR3NyMM2fO+CwjXXl5Oerq6pCeno6oqChERUWhuLgYb731FqKiolzbids5OP3798ewYcN0zw0dOhQ1NTUA+H0OpX/6p3/CCy+8gIceegi33HILsrOzMW/ePBQUFADgtjZCqLap3W7HyZMnPT7/u+++C3q7M8EJgZiYGKSnp6OoqEj3fFFREcaMGROmWpmPpmmYPXs2tmzZgv/5n/9Bamqq7vXU1FTY7Xbddm5ubkZxcbFrO6enpyM6OlpXpra2FgcPHuS++Juf//znOHDgACorK11/GRkZeOSRR1BZWYnBgwdzO4fA3/3d33ksc/DFF19g4MCBAPh9DqXz588jIkLfnEVGRrouE+e2Dr1QbdPMzEw4nU7s3r3bVWbXrl1wOp3Bb/egpiiTS/tl4mvXrtUOHz6s5ebmaj169NC++eabcFfNNJ555hnNZrNpO3bs0Gpra11/58+fd5V5/fXXNZvNpm3ZskU7cOCA9vDDD3u9LPG6667TPvroI23v3r3a3//934u+1LMz3K+i0jRu51DYvXu3FhUVpb322mvaX//6V23jxo1afHy8tmHDBlcZbufQeOyxx7Rrr73WdZn4li1btMTERO355593leG2DlxjY6NWUVGhVVRUaAC0ZcuWaRUVFa7lT0K1TSdOnKjdeuutWmlpqVZaWqrdcsstvEy8u3n77be1gQMHajExMdrIkSNdlzdT5wDw+rd+/XpXmba2Nu2ll17S7Ha7ZrVatTvuuEM7cOCA7nMuXLigzZ49W+vTp48WFxenTZkyRaupqeniaMylY4LD7Rwaf/7zn7W0tDTNarVqQ4YM0d555x3d69zOodHQ0KDNnTtXGzBggBYbG6sNHjxYW7RokdbU1OQqw20duE8++cTrMfmxxx7TNC102/TUqVPaI488oiUkJGgJCQnaI488op05cybo+ls0TdOC6wMiIiIi6l44B4eIiIiUwwSHiIiIlMMEh4iIiJTDBIeIiIiUwwSHiIiIlMMEh4iIiJTDBIeIiIiUwwSHiIiIlMMEh4i6lfHjxyM3N/eq3//NN9/AYrGgsrIyZHUiIvOJCncFiIjcbdmyBdHR0eGuBhGZHBMcIupW+vTpE+4qEJECOERFRN2K+xDVoEGDsHjxYjzxxBNISEjAgAED8M477+jK7969GyNGjEBsbCwyMjJQUVHh8ZmHDx/GpEmT0LNnTyQlJSE7Oxv19fUAgB07diAmJgY7d+50lX/zzTeRmJiI2tpa4wIlIkMxwSGibu3NN990JS45OTl45pln8PnnnwMAzp07hylTpuDmm29GeXk58vPz8dxzz+neX1tbi3HjxmH48OHYs2cPCgsLcfLkSTz44IMAfkiosrOz4XQ6sW/fPixatAhr1qxB//79uzxeIgoNDlERUbc2adIk5OTkAADmz5+Pf/mXf8GOHTswZMgQbNy4Ea2trVi3bh3i4+Px4x//GN9++y2eeeYZ1/tXrVqFkSNHYvHixa7n1q1bh5SUFHzxxRe46aab8Oqrr+Kjjz7CP/zDP+DQoUPIzs7Gr371qy6PlYhChwkOEXVrt956q+uxxWKB3W5HXV0dAKCqqgq33XYb4uPjXWUyMzN17y8vL8cnn3yCnj17enz2V199hZtuugkxMTHYsGEDbr31VgwcOBDLly83Jhgi6jJMcIioW+t4RZXFYkFbWxsAQNO0K76/ra0Nv/jFL7BkyRKP19yHoEpKSgAAp0+fxunTp9GjR49gqk1EYcY5OERkWsOGDcO+fftw4cIF13NlZWW6MiNHjsShQ4cwaNAg3HDDDbq/9iTmq6++wrx587BmzRqMHj0a06dPdyVRRGROTHCIyLSmTZuGiIgIPPnkkzh8+DC2bduGpUuX6sr85je/wenTp/Hwww9j9+7d+Prrr7F9+3Y88cQTaG1tRWtrK7Kzs5GVlYXHH38c69evx8GDB/Hmm2+GKSoiCgUmOERkWj179sSf//xnHD58GCNGjMCiRYs8hqKSk5Pxf//3f2htbcWECROQlpaGuXPnwmazISIiAq+99hq++eYb1+Xndrsd7777Ll588UWuhkxkYhatM4PYRERERCbCHhwiIiJSDhMcIiIiUg4THCIiIlIOExwiIiJSDhMcIiIiUg4THCIiIlIOExwiIiJSDhMcIiIiUg4THCIiIlIOExwiIiJSDhMcIiIiUs7/B8HUXzKnTqUgAAAAAElFTkSuQmCC",
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"(\n",
" particle_df_from_direct_selection[\"temperature\"]\n",
" - particle_df_from_parcels[\"temperature\"]\n",
").plot();\n",
"\n",
"plt.show() # why is this necessary???"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "parcels",
"language": "python",
"name": "parcels"
},
"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.11.4"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment