1 files changed, 429 insertions, 0 deletions
diff --git a/notebooks/plot_barycenter_lp_vs_entropic.ipynb b/notebooks/plot_barycenter_lp_vs_entropic.ipynb
new file mode 100644
index 0000000..e188875
--- /dev/null
+++ b/notebooks/plot_barycenter_lp_vs_entropic.ipynb
@@ -0,0 +1,429 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [],
+   "source": [
+    "%matplotlib inline"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "\n",
+    "# 1D Wasserstein barycenter comparison between exact LP and entropic regularization\n",
+    "\n",
+    "\n",
+    "This example illustrates the computation of regularized Wasserstein Barycenter\n",
+    "as proposed in [3] and exact LP barycenters using standard LP solver.\n",
+    "\n",
+    "It reproduces approximately Figure 3.1 and 3.2 from the following paper:\n",
+    "Cuturi, M., & Peyré, G. (2016). A smoothed dual approach for variational\n",
+    "Wasserstein problems. SIAM Journal on Imaging Sciences, 9(1), 320-343.\n",
+    "\n",
+    "[3] Benamou, J. D., Carlier, G., Cuturi, M., Nenna, L., & Peyré, G. (2015).\n",
+    "Iterative Bregman projections for regularized transportation problems\n",
+    "SIAM Journal on Scientific Computing, 37(2), A1111-A1138.\n",
+    "\n",
+    "\n",
+    "\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "Elapsed time : 0.010513782501220703 s\n",
+      "Primal Feasibility  Dual Feasibility    Duality Gap         Step             Path Parameter      Objective          \n",
+      "1.0                 1.0                 1.0                 -                1.0                 1700.336700337      \n",
+      "0.006776453137632   0.006776453137633   0.006776453137633   0.9932238647293  0.006776453137633   125.6700527543      \n",
+      "0.004018712867874   0.004018712867874   0.004018712867874   0.4301142633     0.004018712867874   12.26594150093      \n",
+      "0.001172775061627   0.001172775061627   0.001172775061627   0.7599932455029  0.001172775061627   0.3378536968897     \n",
+      "0.0004375137005385  0.0004375137005385  0.0004375137005385  0.6422331807989  0.0004375137005385  0.1468420566358     \n",
+      "0.000232669046734   0.0002326690467341  0.000232669046734   0.5016999460893  0.000232669046734   0.09381703231432    \n",
+      "7.430121674303e-05  7.430121674303e-05  7.430121674303e-05  0.7035962305812  7.430121674303e-05  0.0577787025717     \n",
+      "5.321227838876e-05  5.321227838875e-05  5.321227838876e-05  0.308784186441   5.321227838876e-05  0.05266249477203    \n",
+      "1.990900379199e-05  1.990900379196e-05  1.990900379199e-05  0.6520472013244  1.990900379199e-05  0.04526054405519    \n",
+      "6.305442046799e-06  6.30544204682e-06   6.3054420468e-06    0.7073953304075  6.305442046798e-06  0.04237597591383    \n",
+      "2.290148391577e-06  2.290148391582e-06  2.290148391578e-06  0.6941812711492  2.29014839159e-06   0.041522849321      \n",
+      "1.182864875387e-06  1.182864875406e-06  1.182864875427e-06  0.508455204675   1.182864875445e-06  0.04129461872827    \n",
+      "3.626786381529e-07  3.626786382468e-07  3.626786382923e-07  0.7101651572101  3.62678638267e-07   0.04113032448923    \n",
+      "1.539754244902e-07  1.539754249276e-07  1.539754249356e-07  0.6279322066282  1.539754253892e-07  0.04108867636379    \n",
+      "5.193221323143e-08  5.193221463044e-08  5.193221462729e-08  0.6843453436759  5.193221708199e-08  0.04106859618414    \n",
+      "1.888205054507e-08  1.888204779723e-08  1.88820477688e-08   0.6673444085651  1.888205650952e-08  0.041062141752      \n",
+      "5.676855206925e-09  5.676854518888e-09  5.676854517651e-09  0.7281705804232  5.676885442702e-09  0.04105958648713    \n",
+      "3.501157668218e-09  3.501150243546e-09  3.501150216347e-09  0.414020345194   3.501164437194e-09  0.04105916265261    \n",
+      "1.110594251499e-09  1.110590786827e-09  1.11059083379e-09   0.6998954759911  1.110636623476e-09  0.04105870073485    \n",
+      "5.770971626386e-10  5.772456113791e-10  5.772456200156e-10  0.4999769658132  5.77013379477e-10   0.04105859769135    \n",
+      "1.535218204536e-10  1.536993317032e-10  1.536992771966e-10  0.7516471627141  1.536205005991e-10  0.04105851679958    \n",
+      "6.724209350756e-11  6.739211232927e-11  6.739210470901e-11  0.5944802416166  6.735465384341e-11  0.04105850033766    \n",
+      "1.743382199199e-11  1.736445896691e-11  1.736448490761e-11  0.7573407808104  1.734254328931e-11  0.04105849088824    \n",
+      "Optimization terminated successfully.\n",
+      "Elapsed time : 2.89129376411438 s\n",
+      "Elapsed time : 0.014848947525024414 s\n",
+      "Primal Feasibility  Dual Feasibility    Duality Gap         Step             Path Parameter      Objective          \n",
+      "1.0                 1.0                 1.0                 -                1.0                 1700.336700337      \n",
+      "0.006776466288966   0.006776466288966   0.006776466288966   0.9932238515788  0.006776466288966   125.6649255808      \n",
+      "0.004036918865495   0.004036918865495   0.004036918865495   0.4272973099316  0.004036918865495   12.3471617011       \n",
+      "0.00121923268707    0.00121923268707    0.00121923268707    0.749698685599   0.00121923268707    0.3243835647408     \n",
+      "0.0003837422984432  0.0003837422984432  0.0003837422984432  0.6926882608284  0.0003837422984432  0.1361719397493     \n",
+      "0.0001070128410183  0.0001070128410183  0.0001070128410183  0.7643889137854  0.0001070128410183  0.07581952832518    \n",
+      "0.0001001275033711  0.0001001275033711  0.0001001275033711  0.07058704837812 0.0001001275033712  0.0734739493635     \n",
+      "4.550897507844e-05  4.550897507841e-05  4.550897507844e-05  0.5761172484828  4.550897507845e-05  0.05555077655047    \n",
+      "8.557124125522e-06  8.5571241255e-06    8.557124125522e-06  0.8535925441152  8.557124125522e-06  0.04439814660221    \n",
+      "3.611995628407e-06  3.61199562841e-06   3.611995628414e-06  0.6002277331554  3.611995628415e-06  0.04283007762152    \n",
+      "7.590393750365e-07  7.590393750491e-07  7.590393750378e-07  0.8221486533416  7.590393750381e-07  0.04192322976248    \n",
+      "8.299929287441e-08  8.299929286079e-08  8.299929287532e-08  0.9017467938799  8.29992928758e-08   0.04170825633295    \n",
+      "3.117560203449e-10  3.117560130137e-10  3.11756019954e-10   0.997039969226   3.11756019952e-10   0.04168179329766    \n",
+      "1.559749653711e-14  1.558073160926e-14  1.559756940692e-14  0.9999499686183  1.559750643989e-14  0.04168169240444    \n",
+      "Optimization terminated successfully.\n",
+      "Elapsed time : 2.7255496978759766 s\n",
+      "Elapsed time : 0.002989530563354492 s\n",
+      "Primal Feasibility  Dual Feasibility    Duality Gap         Step             Path Parameter      Objective          \n",
+      "1.0                 1.0                 1.0                 -                1.0                 1700.336700337      \n",
+      "0.006774675520727   0.006774675520727   0.006774675520727   0.9932256422636  0.006774675520727   125.6956034743      \n",
+      "0.002048208707562   0.002048208707562   0.002048208707562   0.7343095368143  0.002048208707562   5.213991622123      \n",
+      "0.000269736547478   0.0002697365474781  0.0002697365474781  0.8839403501193  0.000269736547478   0.505938390389      \n",
+      "6.832109993943e-05  6.832109993944e-05  6.832109993944e-05  0.7601171075965  6.832109993943e-05  0.2339657807272     \n",
+      "2.437682932219e-05  2.43768293222e-05   2.437682932219e-05  0.6663448297475  2.437682932219e-05  0.1471256246325     \n",
+      "1.13498321631e-05   1.134983216308e-05  1.13498321631e-05   0.5553643816404  1.13498321631e-05   0.1181584941171     \n",
+      "3.342312725885e-06  3.342312725884e-06  3.342312725885e-06  0.7238133571615  3.342312725885e-06  0.1006387519747     \n",
+      "7.078561231603e-07  7.078561231509e-07  7.078561231604e-07  0.8033142552512  7.078561231603e-07  0.09474734646269    \n",
+      "1.966870956916e-07  1.966870954537e-07  1.966870954468e-07  0.752547917788   1.966870954633e-07  0.09354342735766    \n",
+      "4.19989524849e-10   4.199895164852e-10  4.199895238758e-10  0.9984019849375  4.19989523951e-10   0.09310367785861    \n",
+      "2.101015938666e-14  2.100625691113e-14  2.101023853438e-14  0.999949974425   2.101023691864e-14  0.09310274466458    \n",
+      "Optimization terminated successfully.\n",
+      "Elapsed time : 2.594216823577881 s\n"
+     ]
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcUAAADQCAYAAAB2rXoYAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3Xt8HPV57/HPszdJlmQJW7IxvgMG4hBSiEPIpbmQpDVJCj1NSqDlhKRQTl9tmhtpD7k0TZP2tGnanCav0gsnzUkhF0ICJSYhJUkhp6QJBAMxAXMzNmD5fpNsXVa7O/ucP2ZWWkkrS8arnZH1fb9e+9rdmfHsj2G0zz7P7/ebMXdHREREIBV3A0RERJJCQVFERCSioCgiIhJRUBQREYkoKIqIiEQUFEVERCIKiiIvkJn9k5n9SZ32tcLM+s0sHb3/kZldXY99R/v7npldWa/9iZyoMnE3QCSpzOxZYDFQAgJgM3AjcIO7l939945hP1e7+w8n28bdnwfajrfN0ed9Ejjd3a+o2v9F9di3yIlOmaLI0f2au7cDK4G/Av4n8C/1/AAz049TkYRQUBSZBnfvc/cNwDuBK83sbDP7spn9OYCZdZnZd8ys18wOmtm9ZpYys5uAFcAdUXn0j81slZm5mV1lZs8Dd1ctqw6Qp5nZz8zssJl928wWRJ/1ejPrqW6fmT1rZm8ys/XAR4F3Rp+3KVo/Uo6N2vVxM3vOzPaa2Y1m1hGtq7TjSjN73sz2m9nHqj7nfDPbGLVpj5l9bqaOuUgcFBRFjoG7/wzoAX553Kpro+XdhCXXj4ab+38HnifMONvc/a+r/s3rgBcBvzrJx70L+B1gCWEJ9wvTaN+/A/8L+Eb0eS+tsdm7o8cbgFMJy7Z/P26b1wBnAm8EPmFmL4qWfx74vLvPB04DbpmqTSKziYKiyLHbCSwYt6xIGLxWunvR3e/1qS8s/El3H3D3oUnW3+Tuj7r7APAnwKWVgTjH6beBz7n7VnfvBz4CXDYuS/0zdx9y903AJqASXIvA6WbW5e797n5fHdojkhgKiiLHbilwcNyyzwJbgO+b2VYzu24a+9l+DOufA7JA17RbOblTov1V7ztDmOFW7K56PcjoIKCrgDOAJ8zsATN7Wx3aI5IYCooix8DMXk4YFH9cvdzdj7j7te5+KnAx8CEze2Nl9SS7myqTXF71egVhlrYfGADmVbUpTVi2ne5+dxIOHKredwnYM8W/w92fdvfLgUXAZ4BvmVnrVP9OZLZQUBSZBjObH2VFNwNfcfdfjFv/NjM73cwM6COcwlGOVu8h7Ls7VleY2Vozmwd8CviWuwfAU0Czmb3VzLLAx4Gmqn+3B1hlZpP9fX8d+KCZrTazNkb7IEtTNcjMrjCzbncvA73R4vLR/o3IbKKgKHJ0d5jZEcJS5seAzwHvqbHdGuCHQD/wU+Af3P2eaN1fAh+PRqZ++Bg++ybgy4SlzGbgfRCOhAV+H/gisIMwc6wejfrN6PmAmT1UY79fivb9n8A2IA/84TTbtB54zMz6CQfdXHaUPlGRWcd0k2EREZGQMkUREZGIgqKIiEhEQVFERCSioCgiIhKJ7ULEXV1dvmrVqrg+XkRE5pAHH3xwv7t3T7VdbEFx1apVbNy4Ma6PFxGROcTMnpt6K5VPRURERigoioiIRBQURUREIlMGRTP7UnQj0kcnWW9m9gUz22Jmj5jZefVvpoiIyMybTqb4ZcLrHU7mIsLrPq4BrgH+8fibJSIi0nhTBkV3/08m3juu2iXAjR66D+g0syX1aqCIyIno/tu/yR2f+8u4myHj1KNPcSljb4baEy2bwMyuMbONZrZx3759dfhoEZHZae+2Z9j1zFNxN0PGaehAG3e/wd3Xufu67u4p51CKiJywglKRoFiMuxkyTj2C4g7G3iF8WbRMREQmEZRKlEtT3tdZGqweQXED8K5oFOoFQJ+776rDfkVETlhBsUippEwxaaa8zJuZfR14PdBlZj3AnwJZAHf/J+BO4C3AFmCQ2nclFxGRKsoUk2nKoOjul0+x3oE/qFuLRETmgKBYpBwEeLmMpXQdlaTQ/wkRkRgEUelUJdRkUVAUEYlBEJVOVUJNFgVFEZEYVKZjaFpGsigoiojEoFI+DZQpJoqCoohIDCrBUJlisigoiojEYKR8qoE2iaKgKCISg7LKp4mkoCgi0mDurvJpQikoiog0WHV2qKCYLAqKIiINVq7qR1T5NFkUFEVEGqxUrA6KyhSTREFRRKTBqq9io6CYLAqKIiINVh0I1aeYLAqKIiINNrZ8qj7FJFFQFBFpMJVPk0tBUUSkwapLpiqfJouCoohIg5U0JSOxFBRFRBqsrMn7iaWgKCLSYIHmKSaWgqKISIONucybyqeJMq2gaGbrzexJM9tiZtfVWL/CzO4xs4fN7BEze0v9myoicmLQPMXkmjIomlkauB64CFgLXG5ma8dt9nHgFnc/F7gM+Id6N1RE5EQRaJ5iYk0nUzwf2OLuW929ANwMXDJuGwfmR687gJ31a6KIyIkl0DzFxJpOUFwKbK963xMtq/ZJ4Aoz6wHuBP6w1o7M7Boz22hmG/ft2/cCmisiMvtVMsVMrknl04Sp10Cby4Evu/sy4C3ATWY2Yd/ufoO7r3P3dd3d3XX6aBGR2aWSKeZaWlQ+TZjpBMUdwPKq98uiZdWuAm4BcPefAs1AVz0aKCJyoqmUTHPNLcoUE2Y6QfEBYI2ZrTazHOFAmg3jtnkeeCOAmb2IMCiqPioiUkMlEGZbWigHyhSTZMqg6O4l4L3AXcDjhKNMHzOzT5nZxdFm1wK/a2abgK8D73Z3n6lGi4jMZkGphFmKTC435o4ZEr/MdDZy9zsJB9BUL/tE1evNwKvr2zQRkRNTUCqSzmbJZLIqnyaMrmgjItJgQbFIOpshnc2OuQ6qxE9BUUSkwYJSkXQmSyqTGXPHDImfgqKISIMFxRLpjMqnSaSgKCLSYGGfYoZUJqPyacIoKIqINFilfJrOZlU+TRgFRRGRBguKYVDMaKBN4igoiog0WFAqkc6E5VP1KSaLgqKISINV5immNdAmcRQURUQaLCiVwsn72awuCJ4wCooiIg0W9imG5VP3MuUgiLtJElFQFBFpsMpAm3QmG77XCNTEUFAUEWmwclAa6VOEcDK/JIOCoohIg5Wi8mk6q0wxaRQURUQabKR8mg1vVKSgmBwKiiIiDVYulUhlMlXlUwXFpFBQFBFpsFLVPEVA0zISREFRRKTByqUSmWyWdCYqnypTTAwFRRGRBvJyOC8xldZAmyRSUBQRaaDKXTFUPk0mBUURkQaq3BUjncmofJpA0wqKZrbezJ40sy1mdt0k21xqZpvN7DEz+1p9mykicmKoBMB0NqvyaQJlptrAzNLA9cCbgR7gATPb4O6bq7ZZA3wEeLW7HzKzRTPVYBGR2axUCYqZqoE2Kp8mxnQyxfOBLe6+1d0LwM3AJeO2+V3genc/BODue+vbTBGRE8OY8mlW8xSTZjpBcSmwvep9T7Ss2hnAGWb2X2Z2n5mtr7UjM7vGzDaa2cZ9+/a9sBaLiMxigQbaJFq9BtpkgDXA64HLgf9jZp3jN3L3G9x9nbuv6+7urtNHi4jMHpUAmNY8xUSaTlDcASyver8sWlatB9jg7kV33wY8RRgkRUSkyshAG10QPJGmExQfANaY2WozywGXARvGbXM7YZaImXURllO31rGdIiInhEADbRJtyqDo7iXgvcBdwOPALe7+mJl9yswujja7CzhgZpuBe4A/cvcDM9VoEZHZKtBAm0SbckoGgLvfCdw5btknql478KHoISIik6geaJNKK1NMGl3RRkSkgUaCYiaLmZHOZNSnmCAKiiIiDVTdpwhhxqjyaXIoKIqINFD1lAwIg6PKp8mhoCgi0kCj5dPMyLMyxeRQUBQRaaDqC4JXnsvqU0wMBUURkQYa36eYymQpqXyaGAqKIiINVD1PESCj8mmiKCiKiDTQhD5FlU8TRUFxDnF3HvjuNvoP5eNuisicFZRKpNIZLBV+/aYyWc1TTBAFxTmk/9AwP7tjG888rNt2icQlKBZHskSATDajKRkJoqA4hxSGwj+8Yj6IuSUic1dQKo6MPIUoU1SfYmIoKM4hxeEgetavUpG4BKXSmEwxraCYKAqKc0ghX4qelSmKxCUojs0U01ld0SZJFBTnkMJQGAwrwVFEGi/MFKuCoi4InigKinNIpWyqPkWR+IwfaKPyabIoKM4hlbKpyqci8Rk/0CbMFFW9SQoFxTmkkiEWVT4Vic2E8mlW8xSTREFxDhkpnw4rUxSJy4TyqQbaJIqC4hwyMtBmSH+AInGpWT4tFnH3GFslFQqKc0ghyhQLyhRFYhOUSqTGDbQBKAf6sZoECopzyEif4nCgX6UiMQmKRTLj+hQBlVATYlpB0czWm9mTZrbFzK47ynZvNzM3s3X1a6LUy8ioU1e/okhcyqXShPIpoGkZCTFlUDSzNHA9cBGwFrjczNbW2K4deD9wf70bKfVRPepUcxVF4lEqFWuWT5UpJsN0MsXzgS3uvtXdC8DNwCU1tvs08BlA9yVKqOr5ibqqjUg8gmKRjDLFxJpOUFwKbK963xMtG2Fm5wHL3f27R9uRmV1jZhvNbOO+fbp9UaMV8yVa2sM/RpVPReJRHj/QZqRPUUExCY57oI2ZpYDPAddOta273+Du69x9XXd39/F+tByjwnBAa2dT+FrlU5FYhPMUNdAmqaYTFHcAy6veL4uWVbQDZwM/MrNngQuADRpskyxedor5qqCouYoisQhKJZVPE2w6QfEBYI2ZrTazHHAZsKGy0t373L3L3Ve5+yrgPuBid984Iy2WF6RYCDPD1o4wKKp8KtJ45SDAvayBNgk2ZVB09xLwXuAu4HHgFnd/zMw+ZWYXz3QDpT4qo00rmaKufyrSeJVscPyto6rXSbwyU28C7n4ncOe4ZZ+YZNvXH3+zpN4qo03b1KcoEptKNpjJ1upTVFBMAl3RZo6olEtb2rOYqXwqEodK4EuNyRRVPk0SBcU5ojKwJtecIduU1kAbkRiMlE+z1X2KKp8miYLiHFEpl2ab02SbM7oouEgMKpli7SkZCopJoKA4R1TKpbnmDLnmtAbaiMSg9kAbBcUkUVCcIypBsJIp6tqnIo1X6TdM1xpoU9QP1SRQUJwjKuXTSqaoa5+KNN5o+bTG/RSVKSaCguIcUciXMINMLhUOtFGmKNJwR5unWNJAm0RQUJwjivmAbHMGMyPXovKpSBxGyqe6IHhiKSjOEYXhgFxzGoBcU5rCsMqnIo02Uj6t6lO0VArMKGueYiIoKM4RxXyJbFMYFCsDbdw95laJzC21BtqYGZlMVuXThFBQnCMKUfkUwhGo5cAJSuWYWyUyt4z2KY69wmY6m1WmmBAKinNEMV8aLZ9GwVH9iiKNVWugDUAqk1GfYkIoKM4RhXwwEgwrwVEjUEUaq1b5tPJe5dNkUFCcI8LRp5U+xfC5qME2Ig1Va54iQCaj8mlSKCjOEYV8iVxTZfRp+AdZGFKmKNJIoxcEr1E+VaaYCAqKM8Ddue2zD7Llwb1xNwUI21McN9AGSMxVbYaOFPjqn97HgR39cTdFZEbVmqcIYZAMgmT8Pc51CoozIN9fZNczfex86lDcTQEgKJUpl51cy7iBNgm5U8b+nn569wyye2tf3E0RmVEj91NMjwuKyhQTQ0FxBgz0DUfPhZhbEqqMMs02jc0UkzL6NGnHS2SmBKUS6WwWMxuzPJ3JKigmhILiDOg/FH7J9/cOx9ySUKVMOjolI1nl04HoOA0k5HiJzJSgWJxQOoWofKqBNomgoDgDBqOMZ7AvGV/y1TcYBkaubJOUKRkDveHxGkjI8RKZKWFQzE5YntY8xcSYVlA0s/Vm9qSZbTGz62qs/5CZbTazR8zsP8xsZf2bOntUMsSBvgLlcvyXUitW3TYKIJVOkcmlEnOjYWWKMldUyqfjpTPKFJNiyqBoZmngeuAiYC1wuZmtHbfZw8A6dz8H+Bbw1/Vu6GxSyXi87Awdib+frFB1g+GKbHOGQkIG2oz0KSooygkuKE1SPtVAm8SYTqZ4PrDF3be6ewG4GbikegN3v8fdB6O39wHL6tvM2aX6y30wAYNHKqNMK/MTw9fp5Ay0iY7X0JEiQaDrscqJa9LyaTar8mlCTCcoLgW2V73viZZN5irge7VWmNk1ZrbRzDbu27dv+q2cZQZ6h2lpD0/8JAy2KQzVyhTTiRho42VnoK8wcryS8CNCZKYEpZIG2iRcXQfamNkVwDrgs7XWu/sN7r7O3dd1d3fX86MTZaB3mEUr54+8jttIpthSlSk2J+NGw4NHCnjZE3W8RGZKUCpO0qeo8mlSTCco7gCWV71fFi0bw8zeBHwMuNjd5+w3WxCUGTpSpGt5G2bJ+JIfGX3aNJop5hKSKVYyw0Ur24FkHC+RmVKeNChqnmJSTCcoPgCsMbPVZpYDLgM2VG9gZucC/0wYEJNxbbOYVL7k2xc00zI/l4hpBsV8iUwuRSo1OmE4m5BMsRIEF62KMsUEHC+RmVIqTlI+zWRUPk2IKYOiu5eA9wJ3AY8Dt7j7Y2b2KTO7ONrss0Ab8E0z+7mZbZhkdye8ypd8a2cTrR1Nich8qm8wXJGUPsVKn+vCpW2k0paI4yUyU4420KYclPCyBprFbeJPlhrc/U7gznHLPlH1+k11btesVcl0WjubaO1s4siBfMwtim4wXFU6heSMPh3oGwaD1o5c9CNCA23kxFUuFUnVnLwfLguCgExK11SJk45+nY1kih1hUExCObAwHIwZZAPhoJtSsUw55ikQA73DzGvPkUqnaO1MRrlZZKaUJp28H/59ql8xfgqKdTbQWyCVNlrasrR15sj3FwmK8QaeYj4YM8gGRgfdxH2njIHeAq2dTUCYXat8KieycqlIZpIpGYDmKiaAgmKdDfQOM68jh6WMeR3hl33c2U8hXxq5CHhF5ZJvcV//dKB3eDQoJqQPVmSmBMVJyqcKiomhoFhnA33DtFVlPhD/NIPJBtqE6+IdbDPQNzwmUyzkg9jbJDJTJi+fRkGxqHM/bgqKdTbQO0xrlCFWgmPc9wksHiVTjHOwTVAsk+8v0tqRA0Z/ROiqNnKiKpeKZLK1p2SAMsUkUFCss/HlwMqyOBWPkinGGRSrR+oCI8ExCZfGE5kJk5ZPRzJFBcW4KSjWUSFfopAPRr7km1ozpDOpWINiOShTKpZrZIpR+XQ4vnJN9ZzO6ue4f0SIzAR3j659qj7FJFNQrKNK2a/y5W5mtHbmYs18al3iLXwfDbQZijNTDI9X0vpgRWZCOQh/gNa+ok0lKKpPMW4KinXUPzJHMTeyrLWzicEYR5/Wuhh4+L4yJSMBmWJUZs41Z8g2p2MfrSsyEyqlUc1TTDYFxToaXw6EaJpBjANHRm4wPOGKNvFPyRjoHSadSdHUOhqwdVUbOVFVssCjlU/LyhRjp6BYR+MHjlRe9/cO4+6xtKkykCY3bqBNOpsilbbYB9q0duYwG71QuSbwy4lqJFOsUT5NRctK6lOM3bSufSrTM9A7TLY5PSYAtXY0URoOKOYnXmptJhSDMk/v6efJPYd5YvcRdj9xiDOAD9y2iW23lejPl2jOpmlvznCpG995sIfvMMhZJ7dz5snzOWtJO/ObJ/6SnQnVI3UrWjtz7NrS15DPF2mkyiCaWuXTTFajT5NCQbGOBnoLI/1jFa2do9MMFsxAUHR3tu4f4MdP7+fep/fx02cOMFCIBtekjdc0zQPgjKXzeUlXM61NGYaLZQ7ni7D7IBY43960k6/eH5ZtUgbnLOvktWu6eM2abs5d0Uk2PTMFhYG+Al3L2sYsC8vNYWZdnUGKzHaj5dPJB9qofBo/BcU6CjOf3JhlIyMq+4ZZsKS1Lp/j7jy1p5/vPrKT7/xiF1v3DQCwcuE8fv3cpZy/egEvWjKf1V2tbLl/D3ff+Dgfu+TFzO9qGbOfm5+4nzO7Wvj0772EXX15ntx9hIefP8S9W/bz9/ds4Qt3b6GjJcv6F5/MW89ZwitPW1i3AOnu9PcOs/LFC8csb+1solxy8gNFWtpyk/xrkdnnaANtRsqnyhRjp6BYRwN9wyw5vWPMsnpO4N/ZO8S/PbyDf3t4B1v29pMyuODUhbznVat43RmLWLFw3oR/Uxldmh03TxHCaRnF4QAz45TOFk7pbOENZy3iQ79yJn1DRX76zH6+/9gevvuLXXxj43Y652V560uW8PaXLePc5Z3HlckV8wGl4WBi+bTqeCkozn66P+CoUiEcQFZroE1G8xQTQ0GxTtw9HDgyoXx6fEFxYLjEXY/t5taHevjJMwdwh5evOolP//rZrH/xyXS3Nx3131dGl1ZGm1bLNafJD9Yu13S0ZFl/9hLWn72EfDHg3qf3c8emndz6UA9fvf95Tu1q5TfOW8p/O28ZSztbau7jaEamr4zLrNtOqhyvAl3Ljnm3Uk/uMHQIBvaHz5VHoT96DEJxEMoBlEvhw8uQykAqza5bH6P3v7bF/V+RGAfnNcOapTA0NGFdJVNU+TR+Cop1kh8oUi75hMwn25Qm15I5pmkZ5bJz39YD3PrQDr736C4GCwErFszj/W9cw2+cu6xmRjiZYr5EKm2ksxPLntnmNEcOTn0T5OZsmjevXcyb1y7mSL7I934RBum/+f5T/O0PnuKC1Qt5+8uWcdHZJ9NaI/jWUmukLsC8aI6n5irOsErAO7QNep+Hvh7o2wF92+HILjiyB/r3QHmKzCU7D1JZSKXDh6WgHFDoDej9STOtpxRpWTDFOZZphlwb5FqheT40zY+eO6ClE5o7wn3PcoX9e2DzQwz+4IfwmteNWVcpqap8Gj8FxToZPxG92nSmGbg7j+44zB2P7OQ7m3aysy9Pe1OGi196Cr9x3jJevuqkF1SuLOSDCdMxKnLNmWOep9jenOXSly/n0pcvZ/vBQW57aAe3PdzDh7+5iT+5/VHetHYxv3bOEl53ZjdNmcm/yCY7Xkm5XuwJwT0McAeegYNbRx+HtsGh52D48Njtc20wfynMPwW6zoC2xeGjtRtaTooenWHQys2DTAtMcpf4/R/9GJb7Lqd8414yCzqhMBA+hg+HwXjwIAwdhP69YfDt3wNHdo8GZQ+gDAwAAxa26aRVsGA1LDi16nEaNLXVbEPS9P/8Qdj8EAM/+AHBB68l3dk5si6tTDExFBTrpDLhvFL+q9bakaO/N89gcXDM8qDsbOo5xD1P7OOHj+/huQODZNPGq0/v4gO/spI3vWgxzdkwsAyVJpZcpmNocJhMU2rCZwNYzinkSzXXTcfCdvjd1y3l6teewsPbD7Fh0y6+/1gPdzyyjfbmDG84axEXnrmIV53eRdu4DLL3QH/YhtZgwuc3t2XoOzjwgtt1Ikqn0jSla5TK3cPy5sFnouD3DBzYAgeiAFgcGN02lR0NLCteBSetDN93roCO5WFGVocRv4WeHvq+/W1O+u3fItPdHS5s6QwfLJ16B0EJ+ndD73bofS4M4Ie2waFn4ekfhAG0WtvJsPB0WBgFyYWnhc8LVkP22Ev7M2VkusXQEAdvvJHu971vZF0qlcZSKfUpJoCCYp1UMpt5HRMHh/Rm9rN722Fe8bXfmXwH3dAefX88CDz4KHz60eNv169uuYr5w1284mu/N2Hdy59/C+fl38wrvvoKqNfshxXQHr28Jw/3bAI2Tdzs1dvezhnpdbzm1ldNWPeO8h/xxOaf88GvXV6nRs1+aUvxgVUX8+7W08LgcHBblPltg8KR0Q1TGehcGQaG1b8cZlOVINGxrCFlyAP/fAOWTrPwqqtf2A7SmbCtHctg5Ssnrh/uD4PkyI+A6IfAE3fC4P6x285fGmWUq8MfACOP1WHm28BpP5WA1/7KV3Hwpq+w4N3vJj1//sj6dDar8mkCTCsomtl64PNAGviiu//VuPVNwI3Ay4ADwDvd/dn6NjXZRvrIqsqBg4USX/r5rfz80M95aeFCmg9fQu9QeNK35FKc1tXGmsVtrO5qoyU7M19WQ9sXQZNx7cuunbCuYO0Ud6T40DkfxnL1v+JO2Z3thwbZsqefp/f2c2AgzKabMym6CmsptxiXrPgfnDy/ibam0RF5+Z0LWTS4qGabTzheHi0rDh+BfF/4GD4MQ72Q74V8Hxubm/jbbbcTHDzEVUfyYXa34FRY8crRUuLC08Ll6cZcfKGW4s6d9N5+Oyf95jvILl40Mx/S1AYnvyR8jJfvm1guPrgVnvweDOwbu222NTxencvDADx/afjoWBpmn+2Lw1JxnQJnZZ7iwt+6nEN3/4iDX/kK3b//+yPr05mMyqcJMGVQNLM0cD3wZqAHeMDMNrj75qrNrgIOufvpZnYZ8BngnTPR4CQpl52DgwX2Hh7mqWd7seYUf/G9J3j2wADP7h9ge/HHNC25hbV+CSlPc/68i3nZud284tQFnHXyfNKpmf+V+s07HqC5I8uvnX3RhHWPHtzB//vJk1x62mU1+0LrbWfvED/bdpD7th4g99x+9gQB37xrNQDd7U2c1t3K6q5WlmfypAaGOGf+JSye30x3WxO5TMKvSOgOxaEwsA0fHg1u1aM2Bw+GmczAvrDk2b83fO/jpi2kMtB+SvhlvfjF0LmCKzqW8tE9P+LveBC/8A+5+qXXxPPfOYUDX/wiAAuvfoFZ4vFq7oCl54WP8Yb7w3Lswcrgou3hc+9z0LMx7OMcL9MSBsd5XWHfamsXzFsI8xZAy4Kx/ayVQUJN7TV/mFTKp/POOovihRdy8F9vZMG7riTdFs5fTmeyuqJNAthU1+Q0s1cCn3T3X43efwTA3f+yapu7om1+amYZYDfQ7UfZ+bp163zjxo3H/R9w00f/mKAw+utqsg/0Gu/cw0fVEirTqsruOE7Zoezh91bZywQe9gWWy2VKVd9l8zgNI01f6idh1pctsD29na5MBxfkXsvTu5Zy5pIDtOQae9I/sbOL+fPynLNi/4R1uw+18mhPN2uX7ieXaewv1M09XZzUNsiSk3bQO1Skb6jIwHCJgeESHiwjm1pBPnhsZPtM2simjEwKMuk0aQvHeKTNSKWMFOHVeFIW/rBPMfqu9DCyAAAHRklEQVQD3wCjHD1XlgWYQ/h/3jEcyuVwO3cgel0uh6/L5ejflKEcYJSwcpmUl7BygHlx8pMv4qkMQaqJcjqHp5sop5sop5sJ0s2U000EmWbK6XkE6aaa5Wx32HbkUQ4O72bJvFW0ZTsnbhQjC8oseuARjqw8hQO/9OK4m3PsPCAd5EmVBkkFw6RLeVJBPlwWFEiVh0eeJ/yQGb8rS+GWxVNp3DK4pRg8EnBkX56l555CUz90332AI2vayS9pxS3F889vJ5vNMr9zQbSX1MjVqX3kMtWjJ3VlmY85V4724/FYf4RbzZeNlGnKcsWff6Yu+zKzB9193ZSfOY19LQW2V73vAV4x2TbuXjKzPmAhMOab2MyuAa4BWLFixTQ+emr7tm7FfeppBfVijB60sb2H23FgftWSU1kMwGZ+BsCjT89482o6cgh27Jh8/aYnG9eWagN90FOjXcZeSjw04eQsRo9kmm7JshA9+l/wpywGyjzF4ak2jkHfQqD/WfjxszG35IWZ3vmVjh5TKUeP6r2m2b/t7Zhl2XlOtKg32tq/Qn5oL/mho/yxzjFm059+Vi8NHWjj7jcAN0CYKdZjnytfdh5BMcpybMzTpD90xmYQYJbCDAwjnapkGnbMv46aWwqkM+EvyBTGyU0LSUdD1o8MpCgG8fzc6mgLqHV1NvewXaXy8bTrhf1bA+a3l8a1y0badfhImmDSfY/N+SovgyirD8qj750o06/ki2ZhdcBSUX5oYCnKpMI5dpXdTvPsjOPuJ45zuHCI8hTZShyClmbK89un3nCOSrW0kmkPv+htuEBu72g/Z1B4LcXBfsCrMtEA3KOqRnnkmZEqB2GVg+qS19hzcnT9mLLY0dXttD6+HWWaZr5bZ8JnTmObHcDyqvfLomW1tumJyqcdhANuZtzb/+ijjfgYERGZA6YzeuEBYI2ZrTazHHAZsGHcNhuAK6PX7wDuPlp/ooiISBJNmSlGfYTvBe4iLKR/yd0fM7NPARvdfQPwL8BNZrYFOEgYOEVERGaVafUpuvudwJ3jln2i6nUe+M36Nk1ERKSxEj75S0REpHEUFEVERCJTTt6fsQ822wc8V6fddTFuTqQclY7X9OlYHRsdr2Oj43Vsjud4rXT37qk2ii0o1pOZbZzOlQokpOM1fTpWx0bH69joeB2bRhwvlU9FREQiCooiIiKREyUo3hB3A2YZHa/p07E6Njpex0bH69jM+PE6IfoURURE6uFEyRRFRESOm4KiiIhIZFYHRTNbb2ZPmtkWM7su7vYkjZktN7N7zGyzmT1mZu+Pli8wsx+Y2dPR80lxtzVJzCxtZg+b2Xei96vN7P7oPPtGdGF8Acys08y+ZWZPmNnjZvZKnV+1mdkHo7/DR83s62bWrHNrLDP7kpntNbNHq5bVPJ8s9IXo2D1iZufVow2zNiiaWRq4HrgIWAtcbmZr421V4pSAa919LXAB8AfRMboO+A93XwP8R/ReRr0feLzq/WeA/+3upwOHgKtiaVUyfR74d3c/C3gp4XHT+TWOmS0F3gesc/ezCW+ucBk6t8b7MrB+3LLJzqeLgDXR4xrgH+vRgFkbFIHzgS3uvtXdC8DNwCUxtylR3H2Xuz8UvT5C+IW1lPA4/Wu02b8Cvx5PC5PHzJYBbwW+GL034ELgW9EmOl4RM+sAXkt4lxzcveDuvej8mkwGaInuOTsP2IXOrTHc/T8J77RUbbLz6RLgRg/dB3Sa2ZLjbcNsDopLge1V73uiZVKDma0CzgXuBxa7+65o1W5gcUzNSqK/A/4YqNz6fCHQ6+6l6L3Os1GrgX3A/43KzV80s1Z0fk3g7juAvwGeJwyGfcCD6NyajsnOpxmJAbM5KMo0mVkbcCvwAXc/XL0uuhm05uUAZvY2YK+7Pxh3W2aJDHAe8I/ufi4wwLhSqc6vUNQPdgnhD4lTgFYmlgllCo04n2ZzUNwBLK96vyxaJlXMLEsYEL/q7rdFi/dUygzR89642pcwrwYuNrNnCcvxFxL2mXVGJS/QeVatB+hx9/uj998iDJI6vyZ6E7DN3fe5exG4jfB807k1tcnOpxmJAbM5KD4ArIlGb+UIO603xNymRIn6w/4FeNzdP1e1agNwZfT6SuDbjW5bErn7R9x9mbuvIjyf7nb33wbuAd4RbabjFXH33cB2MzszWvRGYDM6v2p5HrjAzOZFf5eVY6Vza2qTnU8bgHdFo1AvAPqqyqwv2Ky+oo2ZvYWwDygNfMnd/yLmJiWKmb0GuBf4BaN9ZB8l7Fe8BVhBePuuS919fOf2nGZmrwc+7O5vM7NTCTPHBcDDwBXuPhxn+5LCzH6JcFBSDtgKvIfwx7bOr3HM7M+AdxKOCn8YuJqwD0znVsTMvg68nvAWUXuAPwVup8b5FP24+HvCMvQg8B5333jcbZjNQVFERKSeZnP5VEREpK4UFEVERCIKiiIiIhEFRRERkYiCooiISERBUUREJKKgKCIiEvn/SJQF5bnR4kYAAAAASUVORK5CYII=\n",
+      "text/plain": [
+       "<Figure size 460.8x216 with 1 Axes>"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZUAAAELCAYAAAARNxsIAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzsnXeYVNXdxz9ndrZXtrK0pXeWhQVpNkTFWACNxBo1dk2i0cRoTOKL7bUkYntVEpWYRKNGE40oMQYLFgQEBKSXpS0sW9neppz3jzt3drZO2Tszd3fP53nmmdk7955zZvbO/d5fOb8jpJQoFAqFQmEElnAPQKFQKBS9ByUqCoVCoTAMJSoKhUKhMAwlKgqFQqEwDCUqCoVCoTAMJSoKhUKhMAwlKooehxBimRDitwa1NUQIUSuEiHD9/ZkQ4noj2na1928hxNVGtadQmB1ruAegULRFCHEQyALsgAPYAfwF+KOU0imlvNmPdq6XUq7qbB8p5WEgobtjdvW3BBgppbzSo/3vGdG2QtFTUJaKwqxcIKVMBHKAR4G7gZeN7EAIoW6qFAqDUaKiMDVSyiop5XvAJcDVQoiJQohXhBAPAQgh0oUQ7wshKoUQFUKIL4QQFiHEX4EhwAqXe+uXQoihQggphLhOCHEY+MRjm6fAjBBCrBdCVAsh/iWESHX1dboQotBzfEKIg0KIM4UQ5wD3Ape4+tviet/tTnON6zdCiENCiBIhxF+EEMmu9/RxXC2EOCyEKBNC/Nqjn5OEEBtcYyoWQiwN1neuUHQHJSqKHoGUcj1QCJzS5q2fu7ZnoLnM7tV2lz8EDqNZPAlSysc9jjkNGAfM76S7q4BrgWw0F9wzPozvQ+B/gTdd/U3uYLdrXI+5wHA0t9v/tdnnZGAMMA+4TwgxzrX9aeBpKWUSMAL4u7cxKRThQImKoidxDEhts82GdvHPkVLapJRfSO8F7ZZIKeuklA2dvP9XKeU2KWUd8FvgB3ogv5tcASyVUhZIKWuBXwGXtrGS7pdSNkgptwBbAF2cbMBIIUS6lLJWSrnWgPEoFIajREXRkxgIVLTZ9jtgH/CREKJACHGPD+0c8eP9Q0AkkO7zKDtngKs9z7ataBaWznGP1/W0JBFcB4wGdgkhvhFCnG/AeBQKw1GiougRCCGmo4nKl57bpZQ1UsqfSymHAwuAO4UQ8/S3O2nOmyUz2OP1EDQroQyoA+I8xhSB5nbztd1jaIkHnm3bgWIvxyGl3CulvAzIBB4D3hZCxHs7TqEINUpUFKZGCJHkuit/A3hVSvldm/fPF0KMFEIIoAotBdnpersYLXbhL1cKIcYLIeKAB4C3pZQOYA8QI4Q4TwgRCfwGiPY4rhgYKoTo7Hf1OnCHEGKYECKBlhiM3duAhBBXCiEypJROoNK12dnVMQpFOFCiojArK4QQNWiuqF8DS4EfdbDfKGAVUAt8DTwvpfzU9d4jwG9cmWG/8KPvvwKvoLmiYoDbQMtEA24FXgKOolkuntlgb7mey4UQmzpod7mr7c+BA0Aj8FMfx3QOsF0IUYsWtL+0i5iQQhE2hFqkS6FQKBRGoSwVhUKhUBiGEhWFQqFQGIYSFYVCoVAYhhIVhUKhUBiGEhWFQqFQGIZfVVrT09Pl0KFDgzQURSjYuHFjmZQyw/ue/qHOjZ5PsM4NUOdHb8DX88MvURk6dCgbNmwIfFSKsCOEOOR9L/9R50bPJ1jnBqjzozfg6/mh3F8KhUKhMAxTLVLUbHfy9sZCzsvNJjk2MtzDUSgCQ0r4+jmo9VrSyzdShsBJNxjTlqL3se9jSMyGrPHhHglgMlF58YsCfvef3Ww4VMHSH+SFezgKRWBUH4OPXOtrWWO715bdVYkl9xKISepeW4reyQd3wpDZcOEL4R4JYCJRKalu5NlP9pIcG8k/Nx3lhzNzmDKkX7iHpVD4j9OmPS98HqZc0b221i6DD+8Gp9eak4q+isMGjuZwj8KNaWIqXxeU02hz8scf5mMR8Nnu0nAPSaEIDKdDe7YYsK6X3obepkLRFqcDpHnOD9OIyuYjlcRGRpCf04/RWYlsKaz0fpBCYUakqyK9EYtF6lX0TXTRUJgM6TDVTYdpRGXLkUomDUzGGmFh8qAUthypRFVQVvRI3JaKAT8vZakovOF0tNzImABTiIrN4WT7sWomD04GYPLgFE7U2zhSoZaLUPRAdKvCEEslonWbCkVblKXSnt3Ha2iyO8kdlALgFpfNygWm6ImomIoilDidprrpMIWo7C+tBWBM/0QARmYmIATsL6kN57AUisDQf+AWA5Ir9TZM5N5QmAxlqbTncHk9AENS4wCItkYwIDmWwxX14RyWQhEYTiPdX66fqEopVnSG02Gq88MUonKoop6spGhiIlt+hENS4zhUXhfGUSkUAaIC9YpQ4rSbypI1hagcLq8nJzW+1bactDhlqSh6JipQrwglyv3VnkMVdQxJi2u1bUhaHGW1zdQ2mcesUyh8QgXqFaHC6bJQTHTTEXZRabQ5KK5uIie1tajoloseb1EoegzKUlGECv28MNFNR9hFRXdxtbVUclx/H65QcRVFDyMolop5fOYKE6Gfaya66Qi/qLTJ/NIZ4hYVZakoehiqTIsiVChLpT1F1Y0ADExpXSI8KSaShGgrRVWN4RiWQhE4KvtLESrclop5LNmwi8rxqgasFkFaQnS79/onx3BciYqip6FiKopQoSyV9hRVNpKVFEOERbR7Lzs5hmNKVBQ9DZX9pQgVKvurPUVVjfRPjunwvf5JMRyvUkUlFT0MZakoQoWyVNpzvLqR7E5EJTsllpKaJmwO8/gLFQqvuC0VA2t/meiioTAR+nmhyrRoSCkpqmroXFSSY5ASSmuaQjwyhaIbGOr+0mt/KVFRdIBUgfpWVNbbaLQ56Z8c2+H7ulusSLnAFD0Jt/vLgJ+Xcn8pusKp3F+t0NOFu7JUPPdTKHoEKlCvCBW628tENx1hFZXj1ZoF0lmgPjtJs2CKKpWoKHoQKlCvCBW628tENx1hFZXiai1W0j+pY1FJirUSE2mhpEaJiqIHocq0KEKFKtPSmmLXbPqMxPYTHwGEEGQlxbjFR6HoEagyLYpQ4U4pNs9NR9gtlfSEKCIjOh9GZmK0W3wUih6BiqkoQoWyVFpTUt1IZmLHri+dzKQYSlRKsaInobK/FKFCTX5sTXFNI5lJHbu+dLISYyhRloqiJ6EsFUWoUGVaWlNS3USWF0slKymaumaHWgFS0XNQ2V+KUKEslRbsDidltU1kebFUdEtGxVUUPYaglGkxTyBWYSI8YypShncsLsImKuV1zTilFjPpCt2SKTFTBlh9BdQUh3sUCrMSlDItylJXdICnBWuSUi1hExXd8sjyIiq66JhmrsonD8ETY2DpWHj3x6YyOxUmQQXqFaHC8/pjkmuRAfZ5YOhzTzI7maOiYyr31zcvw+e/g4kXQ3w6rFsG8Wlw1gPhHpnCTDgdmqCI9msE+Y0K1Cu6opWlYo5zJGyiolse3rK/EqOtxEZGhH8CZF05/Pd/YMQZcNEftR+7rR7WPAt5V0DGmPCOT2EepMOYID0oS0XRNZ6xNpPceITR/dWEEJDRwTLCnmiz6qPDP1dlzdPQXAvzH2m5e5y3BCLj4NOHwzo0hclwOoyJp4Aq06LoGs9Ym0luPMImKiXVjaTFR2PtYja9TmZSTHjdX831sOFPMPEiyBzbsj0+DU66EXa8B1WF4RufwlxIp7JUFKFBmi+mEj5RqfGeTqyTlRTmCZA7/gVN1TDt2vbv5V8NSPj2tZAPS2FSDLVU1CJdii5wquwvN8XVjV4zv3QyEzX3lwxXHva3r0LqcMiZ0/69fkNh+Omw+VXT5Ikrwox0GJP5pSMilKWi6BhlqbRQXN3kNfNLJyspmvpwzaqvKYZDX8GkH3SezZN7CVQehmObQjs2hTkx0lIBrS2TXDAUJsNpvuyvsIiK3eGkvK7J68RHHd2iCUsG2O4PAAnjF3S+z+hztJnPO1eEbFgKE2Nk9hcoS0XROVJlfwFQVtuMlPgcU8l0z6oPQ1xl5wrN9ZU5vvN94lJh6ClawF65wBRBsVTM4S9XmAxlqWjomVzeyt7ruCdAhnpWfVMNHPgCxp7nfSLb2POgYj+U7w/N2BTmxekwpu6XjiVClWnxgU92FVPX1wrPmjCmEpbJj/qcE3+yvyAM9b8OfA5OG4w8y/u+I8/Unvf9F9JHBndcvZAtRyo5VFEf7mF0yrC0eCYNSvZtZxWoDzlltU1c+8oGHvv+JC6ZPiTcwwkdqkyLxnE/LZWEaCvxURHu40LG3v9CVAIMmeV939RhkDZSO2bmLcEfWy/jipfWmXp5g9T4KDb91oebC1CB+jDQ0Kx9P/XNfex7UmVaNIqrGomwiE7Xpu+IrOQQT4CUEvZ9DMNOA2uUb8eMPAs2LAdbA0TGBnd8vQgpJbVNds7LzeaOM0eHezjteOzDXaza6UdVahWoDznNDi3mZHP0sdiTslQ0iqoayUyMJsLie8G97OQYiqpCKCoVBVB1GObc5vsxI86AdS/A4bUwYm7wxtbLsDm05IZx/RMZmZkQ5tG0Z1z/RP67oxgpJcKXIpEqUB9ybG5R6WOJMp7ZXya58QhboL5/sm+uL53+SbEUh1JUDqzWnof7IQ45s7UArX6swif0C0KUNawLkXaKPi6708cLlpFlWkBZKj7QbNfOoSZ7HxNfE1oqYfkVF1U10N/HOSo6/ZOjKa5pwuHrD7u7FHwGSQMhbYTvx0QnwKDp2rEKn9EvCJE+1IELB/q4mn29YDkdLeVVjMBiMc0Fw6zY+qr7Sy3SpXG8KgBLJTkWh1NSVhuCDDCnU8v8Gnaa/2tiDD8djm3WVodU+ITZLRVdVHy+YKmYSshptms3m7Y+Z6l4JLeY5MYj5L/imkYbdc0Osv0UlWyXZXM8FC6w4m3QcAKGn+b/scNOBSQcWmP4sHorepDVrJaKLnbNvoqKyv4KOX3WUlGTH1tEwddikjq6ZROSYP2Bz7Xnoaf4f+zAaWCNbWlD4RXdrRRlVlHx1/2lLJWQo4tJc18O1JvkxiP0ouJKC85O9i/lVheVkKQVH/wCUkdA8kD/j7VGwZCZWhsKn9Azdkzr/rJqLlCfM4tU9lfIUZYKprnxCPmvWLc0/HV/pcZFERVh4VhVQzCG1YLDrrmuhgVgpegMOxVKdkBtqXHj6sXYzO7+itAEwveYitHZXxbTXDDMim6h+GxN9hZMWKYl5L/iY5UNCOF9bfq2WCyC/skxFFUG2VIp2qItyDXs1MDb0I9V1opPtMRU/EyKCBH6uHzP/rIbbKlYVe0vL+j/mz5tqfRVUSk80UBWYgzRVv9/dIP6xXK0MsiWysFuxFN0svMgKlGJio+4YyqmdX+pQL3Z6bPuLxOWaQn5r/joiQYG9gushMnAlFgKTwS56OCBLyBjLCRkBt5GhBVyZmltKbziTik2qfsrWk8pVoF609JnA/XKUoHCynoGBSgqg/rFUVzdRJM9SF+evRkOf909K0Vn6ClQvheqj3W/rV6O2WMquqWiAvXmxe3+6nMxlT5epsXucFJU2dgNUdGOC1pc5egGsNUHNj+lLXobKrXYK6Z3f+kpxQ4ff7QqUB9ydMH32UXZW+jrlkpxTRN2p2RgSlxAx+tus8ITQYqrFKzWfsBDT+5+W1mTIDZVa1PRJbrLwqyWSss8FX8sFSPLtKiYijdUTIW+WablqEsMumupBC2ucmA1ZE+G2H7db8ti0dKSD6xWSwx7wWb2yY/ueSoqpmJW3DGVvub+6uuWii4GgYpK/6QYIiwiOBlgTTVQ+I1W78sohp8O1UehfJ9xbfZCmntI7S//Ckqq7K9Q0mdTivt69tfhinqEgAEpgYmKNcJCdnIMh4Ox7OyBL7S5ACPOMK5NvWz+vo+Na7MXYjP9PBVVUNLstCzS1ce8An3dUjlQVseA5FhiIgP/wQ1Lj+dAWZ2Bo3Kx/2OIjNNKrBhF6jBIHa61regUd+l7k1oqUVY/RcXpVNlfIabPxlScDrBEaq9NcuMRclEZnhHfrTaGp8dzoLQOaXScYt/H2kx4q38z/b0y8kw4+CXYQ1Cyv4fSbPJ5Krql4vMCUIZbKn03++urfWVc86f1OL2so2Sz+16m5fnP9vHER7sNGV/YkQ6IcC133tcsFSklB0rrGJ7ePVEZlh5PTZOdUiPXVSnfDycOwIh5xrWpM2KelqZ86Cvj2+4l6BcEs2d/+T5PRZVpMYp1Byr4bHcpNU1df/6WyY/eReWzXaV8sqvEkPGFHacDIiJbXpuAkP2KS2ubqGmyM6ybojI8Q1vD/ECpgS6w3f/WnkfPN65NnWGngjUGdn9ofNu9BJvDSYRFEGExZ0zFf/eXCtQbRb1LTOq8iEqzH+6v2ia71/Z6DJ6Wikms2ZCJii4CuigEii5KBUbGVfZ8CJnjoV+OcW3qRMVpWWB7/q1Sizuh2eE0resLIMIisAi1nko4qGvWLv71zV5Exe57oL6+2U5dcy/5Pp192P2lB9e7a6kMSIklymoxLljfcEIrdT/6HGPa64jR50DlYa0cvqIdzXanaTO/dCIjLCpQHwbqmhytnjtD/984nBKHl/hLXbPDbQH1eKSzxf1lkhuPkInK3pJaYiItDAwwnVgnwiIYnh7PnuIaYwa2a6X2zxh7vjHtdcSYcwEBO94LXh89GJvDado5KjpRVovvJUBUoN4wdAvFm7vK00LxJv51TZql4i343yNwOlqSi/qapbL9WBVj+ydhMcBvPj47iR3Hqg0YFbD9HUgZAgOnGtNeRyRmQc4crS9FO2wmd3+BFqz3L6aiyrQYQa0eU/HirvIU/K7E3+mU1LvaarD1gu9UOvqmpSKlZPuxaiYMSDKkvQkDkympaaKkppuFJRtOQMFnMH4hiCC7XyYsgrLdUKxcYG1ptjtNO0dFJzLComIqYUAXAG8xFU/B76pSsaeQ1Hlps0fQKqZiDhdpSH7JRyoaqGm0M2FAsiHt6eK0vbvWyrZ/gNMGEy82YFReGL9QuzhsfSP4ffUwbA5p2nRinSirJcyl7/umqOhur1ofU4q1153/nzzdaN7iND0C2UcnP24/VgXAxIHGWCrjXaLSbRfYt69B1kStiGSwScjUUpa3vAGOXnCHZCBmz/4CrYRM+GIqEaapQBtq3JaKt0C93beYiqcbrVekFTud2jwmYTHNjUdIfsnfHa0iwiIYnZVoSHtJMZHkpMXxXWFV4I0Ub4djmyDviuC7vnTyroDaYtj7n9D010PoVe4v3QWhLBVDaImp+DZPpe3rtngKSX1vSCuWrvidiVykIfklry0oJ3dQcrdqfrVlWk4q6w9WBJ7BsfYFsMbC5EsNG5NXRs+HpIGwblno+uwBaIF6c6cUa+4vH0RF/2Gr7K9uI2VLUN3r5Ee7k/go7Tvv0lJp5f7qDZaKyyo20Y1H0EWltsnOlsIqZo9IM7Td2SPSqKhrZtfxAFKLa0th698h7zKISzV0XF0SEQkn3aCtBnn8u9D1a3J6REqxr9lf+g9bZX91mya70z3nxFv2l83hJC7aCnQ9SdXTOukdgXpXSSATuUiD/kv+5kAFDqdk9oh0Q9udPVITqTX7y/w/+KuntAD9zFsNHZNPTL0aohJh9eOh79ukaJMfzS0qvru/XBcqi9W4zvto7a9WriofAvU+WSrNnm32AqGWnpaKOc6RoP+SP99bSpTVQn6OAaspepCdHMvw9Hi+2OunqFQfg29egtxLIX2UoWPyibhUmPVj2PkeHN0U+v5NSHMPyP6KtFrcyx53SVDcXxGA7HNlfjytilqvM+olcVG6peJb9pe3jLIegV69wUTWbFB/yXaHk/e3FjF3TIah8RSdsyf058t9ZZT5U7H4w3u059PvNnw8PjPrxxCfASt/YZoTIZz0GPeXT5aK7v4yOFDv2XYfoZVV4UOgPsHl/uo6ptLyHXprs0cgHVrMra8E6r/YV0ZpTRMXThkUlPYvmjoQh1Py3uZjvh2w/V3Y8S849S7oNzQoY/KJmCSY/wgc3QhfPxe+cZiEZrv5U4qjrD6mFOt+baMD9WCai0ao0AXAIrqOqUgpXTEV7+4vXUi8tdlj0OdE9RVL5a9fHyIlLpK5YzOC0v7orEQmDUzmtXWHsHv7wZfuhn/9BAbmw+zbgjIev5h0sVZvbNUSbSnjPozN0YsKSgYrUO/Zdh9Bd1WlJUR3malld0qkhPgo74H62iYHkRGCpNjI3pH9pcdU+oKlsq6gnE92lXDjqcOJthrv+tL58dyR7C+t4+2NhZ3vVLYX/rIQImNh8Z/BGhW08fiMELDoeUgbAa9fBofXhntEYaNXub+CFlPBNBeNUKFbFZmJ0V0G6nWxj3dZKl1ZlPXNduKjrcRHWXvHjPpWlkovzv46UdfM3f/YSnZyDNfOGRaMLtzMn5BFfk4//nflTvaX1rZ+U0rY+ha8eAY4bHDVvyBlcFDH4xcxydqYEjLhzxdoc2f62N0oaKmjPSNQr2IqoUS/6GcmRnfpqtJn0+uB+q7LtDiIj7ISHx3Ri2IqEaaay2T4L3nHsWoue3Etx6oa+b/LpwQlQO+JEIKnLskjMsLCFS+u44u9pci6ctj8N01M/nk9pI+GGz+DrPFBHUtAJA2A61fB8LlaEsFzM2D9i1BdFO6RhYyeUqXYp5TioFoq5rgTDRV6oD4jUXN/yU6y33Sx9y1QbycuKoK4KKvK/goS3UqmX//+yzQ3NdBoc1LZYONweR2HKxqYFm3lqVOGMLayGiq7O0Q9ldLz2anlZDtsYGtgcHMt/x1znJ27d5H91yMIy3EAKmMGs2PC/1CYcyFiH1hEIUKEriqLX4x7ggEp5zBmzx9JWfkLWPkL6mIHUJM4gobY/jRH9cMWGY/DEoMzIgopIpBYQAAIJPoHa/lwQ6edS/qAIKxm6QNHC7ZTtO1zn/Y9x3GYySd2w5atQR5V4EypPEplcykb3ivocr/4usOMAzYcqeIIXbhk/WDY0WrygO8+eJ6maGMmEUcm9GPyGSGsJtGGzf/9G/aGrsssyWPVLLJUclpjFs0Us/5fBztccvpEvY1FliLyTuxnkaWI4i83seF4x/PihheVMCbCQrTdQmWJjQ3v9ey0/mnVhRysaCStyUn89nfZFDXNkHYH5J7BgKFjAjq2W6IycsMSUmlT1DEKkMDX3WnZXwSpsSnMSs2m0DKBtxrPYUXNaD6vzIGNAjb2lHLzWcBvGS2OcKplK5Md+xlWd4QssYVUaogQ/s1T+C5lQPhEZcvHnLTltz7tOy0K2Ot6mJSFwMIowMdr0NNrq/hizRZD+j7HUseyKJi0/XeGtAewL2IEhFFU+n39MDnOrkV3GnB1FLAPzosCNne+71lRwC44MwqoptP/U7tLbs/WFADeL7AzzRLLTIudaZvuMaTNjRFPhEdU7Nf+l3Ik8VFWYiKD6L7Q78D1Z70sQUSUtupZZBxYLFiAIa7HYqDJ7qCm0U6T3YnTKXFK2SPnj9UBdVIi7A3aw2ED6UBIJ+0tuRZGZodHUADGnXElhZPO8GnfCIugf3IMZjQgdSRQVNmI04cTSFqjeTBhgIG9n86Rmsu1/7tBJETHGNZWIERd/U8Kbd4/T3JcJAnRVo5XN+LoIlYSabWQlRhNdaON6oau3VqZSTEIAcVV3VyPyQwIwflJQxD2Rgrrig1rdkzmwICP7ZaoZA4Z253Dg060NYLohODGdEJLQrgH4DOJyakkJoewrlqQEcAAY8vX+Ud6YHeNZiU7x7/Pk+3jd5/kevjCoODMdAgTiYA5PpC5o6MKhUKh6FEoUVEoFAqFYYjO0vQ63FmIUuBQ8IajCAE5UkrD7WR1bvQKgnJugDo/egk+nR9+iYpCoVAoFF2h3F8KhUKhMAwlKgqFQqEwDCUqCoVCoTAMJSoKhUKhMAwlKgqFQqEwDCUqCoVCoTAMJSoKhUKhMAwlKgqFQqEwDCUqCoVCoTAMJSoKhUKhMAwlKgqFQqEwDCUqCoVCoTAMJSoKhUKhMAwlKgqFQqEwDCUqCoVCoTAMJSoKhUKhMAwlKgqFQqEwDCUqCoVCoTAMJSoKhUKhMAwlKgqFQqEwDCUqCoVCoTAMJSoKhUKhMAwlKgqFQqEwDCUqCoVCoTAMJSoKhUKhMAwlKgqFQqEwDCUqCoVCoTAMJSoK0yKEOCiEaBBC1AohTgghPhBCDA73uHxFCLFECPFquMehUIQSJSoKs3OBlDIByAaKgWf9bUAIYTV8VCGgp45b0bdRoqLoEUgpG4G3gfEAQojzhBDfCiGqhRBHhBBL9H2FEEOFEFIIcZ0Q4jDwicvK+alnm0KIrUKIC12vJwgh/iuEqBBCFAsh7nVttwgh7hFC7BdClAsh/i6ESG3Tz9VCiMNCiDIhxK9d750D3Atc4rK0tri2JwshXhZCFAkhjgohHhJCRLjeu0YI8ZUQ4kkhRDmwRAgxUgixWghR5Wr/zaB+0QpFN1GiougRCCHigEuAta5NdcBVQApwHnCLEGJRm8NOA8YB84E/A1d6tDcZGAh8IIRIBFYBHwIDgJHAx65dfwoscrU1ADgBPNemn5OBMcA84D4hxDgp5YfA/wJvSikTpJSTXfu+AthdfUwBzgau92hrBlAAZAEPAw8CHwH9gEEEYKkpFKFEiYrC7LwrhKgEqoCzgN8BSCk/k1J+J6V0Sim3Aq+jXfg9WSKlrJNSNgDvAaOFEKNc7/0Q7YLfDJwPHJdSPiGlbJRS1kgp17n2uxn4tZSyUErZBCwBLm7jmrpfStkgpdwCbAEm0wFCiCzgXOBnrnGVAE8Cl3rsdkxK+ayU0u4atw3IAQa4xvalf1+fQhFalKgozM4iKWUKEAP8BFgthOgvhJghhPhUCFEqhKhCu/intzn2iP7C5T57E7hSCGEBLgP+6np7MLC/k/5zgHeEEJUucdsJONAsCZ3jHq/rgYQu2ooEijza+wOQ2dGYXfwSEMB6IcR2IcS1nbStUJgCJSqKHoGU0iGl/CfaBf1k4G9o1sdgKWUysAzt4tvqsDZ//xm4As1NVS+l/Nq1/QgwvJOujwDfk1KmeDxipJRHfRl2B201AekebSVJKSd0doyU8riU8gYng8emAAAgAElEQVQp5QDgJuB5IcRIH/pWKMKCEhVFj0BoLESLLewEEoEKKWWjEOIk4HJvbbhExAk8QYuVAvA+kC2E+JkQIloIkSiEmOF6bxnwsBAixzWODNc4fKEYGOqyjJBSFqHFR54QQiS5kgBGCCHauu08P/diIcQg158n0ETH6WP/CkXIUaKiMDsrhBC1QDVa4PpqKeV24FbgASFEDXAf8Hcf2/sLMAlwzx+RUtagxWsuQHNl7QXmut5+Gs0i+sjV11q0YLovvOV6LhdCbHK9vgqIAnagicTbaOnSnTEdWOf6Dt4DbpdSFvjYv0IRcoSUbS10haL3IoS4CrhRSnlyuMeiUPRGlKWi6DO40pJvBf4Y7rEoFL0VJSqKPoEQYj5Qihbn+FuYh6NQ9FqU+0uhUCgUhqEsFYVCoVAYhhIVhUKhUBiGX1VQ09PT5dChQ4M0FEUo2LhxY5mUMsPodtW50fMJ1rkB6vzoDfh6fvglKkOHDmXDhg2Bj0oRdoQQh4LRrjo3ej7BOjdAnR+9AV/PD+X+UigUCoVh9HpRWbUKpk+Hmppwj0RhNi6+GP7yl3CPQtEbOO88eOedMHV+3XXw/PNh6rw9vV5U/v1v2LAhjP9whWl5/334/PNwj0LRG/j3v2H9+jB1vmoVrFkTps7b0+uXK925U3t+9VW46qrwjkVhHpqatIeyYANj48aNmVar9SVgIj7cnD7++OPs1H+MvZCVKyEpqeV6E1Jeegmiow3rPCYmhkGDBhEZGRnQ8X1CVISAjz+GsjJIb7vihqJPUl2tPStRCQyr1fpS//79x2VkZJywWCxeZ1Dv2LEjZ9y4caEYWsiREurqICsLBg8OwwCamyEhAUaM6HZTUkrKy8spLCxk2LBhAbXRq91f9fVw6BCccgo4nbB1a7hHpDALSlS6zcSMjIxqXwSlt6MXJQlbcRIpDetcCEFaWhqNjY0Bt9ErRaWqCrZtg927te/64ou17du3h3dcCvPgi6i89RYsWRKS4fRELEpQNHpbpSsh2q515x+9UlTuuw8mTYLHHtf+248f/R6JyTa2bQvzwBSmoTNRWbMGvvpKe/3GG/DCC9rrd9+Fhb4uzaUICQkJ7VdtXrp0KePHjyc3N5d58+Zx6FDHUysaGho47bTTcDgcHDx4kNjYWPLy8tyPvwSQFihly5iOHTvGxfrdbBtOP/10Y+fsSMl3u3ZxzTXXGNdmN+iVMZXdu7XnN98QEF9MXdJmalPWsmnLdLSlzhV9me++g2PHtNdtReXuu8Fuh6+/hvJyqKjQLharVsF772ku1YoKeOUV+PWvtXidN8rLtSBugHFPhR9MmTKFDRs2EBcXxwsvvMAvf/lL3nzzzXb7LV++nIsuuoiIiAgARowYwebNmwPqsyP314ABA3j77bcDai+QAUwaPZrCwkIOHz7MkCFDQtNvJ/RKS+XECTj5FEnm7Qs47ckfseXH6yFzG99tk73OVFX4R329Nm/pwQe1v9uKSnGx9gBNDOx2zaopKWnZ9sYb8NvfwsGD2raKis77s9u15JCbbzb0Yyg6Ye7cucTFxQEwc+ZMCgsLO9zvtddeY6EPpqenNfT222+7rYHi4mIuvPBCJk+ezNSpk9mypXVK78GDB5k4cSKgWUWXXnop48aN48ILL6ShocG930cffcSsWbOYOnUqixcvpra2FoAHHniA6dOnM3HiRG688Ub0avKnn346d999NyeddBKjR4/miy++cKvZBRdcwBtvvOHL1xRUeqWlcvgwTJxziJJ+K1h28j8ZnDyY3IkRbPkmlj37mxgzMjrcQ1SEid27tVRi3ZptbNQu/FbXL6G0FGw27XV5ectzaWnLa93KKSmBo0fhtNNg1y4YNap9f7oYvfdecD6PKbj22sFs2xbX1S45AHFd7tKavDx46qluDevll1/me9/7Xrvtzc3NFBQU4FmLbP/+/eTl5bn/fvbZZznllFM6bfu2227jtNNO45133qGx0cHatbWd3rC+8MILxMXFsXPnTrZu3crUqVMBKCsr46GHHmLVqlXEx8fz2GOPsXTpUu677z5+8pOfcN999wHwwx/+kPfff58LLrgAALvdzvr161m5ciX3338/qx57DKRk2rRpPProo/zyl7/096sylF4nKk1NcPw4ZEVuJD0unQvGaP+IWy8ZzU1/gqUvH+IPj4wO8ygV4UJP5Xc6W7bV1EC/flpmZmWltk13c4H2rItKWZkmJKAJRkGB1tbu3ZqofP01zJgBFpcPoKhIe87MDO7nUrTm1VdfZcOGDaxevbrde2VlZaSkpLTa5q/765NPPnHHXSyWCBISkjsVlc8//5zbbrsNgNzcXHJzcwFYu3YtO3bsYM6cOYAmdrNmzQLg008/5fHHH6e+vp6KigomTJjgFpWLLroIgPz8fA4ePOi2VDIzMzmm3/GEkV4nKvoP/pD8krOHnYHVon3EK+ZO5+aB6/nnG0O4+Aw480zf/OGK3sWOHe23VVdrolJW1rLt4EHtBgU060S3OMrKWlsq+/Zpr48dg3XrYPbs1kH9PiEqy5cf8bbLoR078sePHx+K0bBq1SoefvhhVq9eTXR0e69EbGyszymznplQ3Umz7QgpJWeddRavv/56q+2NjY3ceuutbNiwgcGDB7NkyZJWfeufKSIiArvdrjdGY2MjsbGxho4xEHpdTOXwYe25Mnor84bNc2+Pj4pn+OlrKDvYn7PP1kq3KPoeHYmKHlfRrRHQ3Fk6paUtrrCuROXjj7XXnhOb9X0zglJQXtGWb7/9lptuuon33nuPzE6UvF+/fjgcDp9EIisri507d+J0OnnHo9bTvHnzeMGVGmi3O6itrerUUjn11FP529+0Fay3bdvGVteEuZkzZ/LVV1+xz3US1dXVsWfPHve40tPTqa2t9S3gLyV79uxxx3HCSa8TlSP6PVPy4VaiAnDZ1TVw8qOt91P0KTqqZNGRqHjut3dvi7usrajs36+9PnaspY5YQUGLW0y3VFQlB+Opr69n0KBB7sfSpUu56667qK2tZfHixeTl5bFgwYIOjz377LP58ssv3X/rMRX98cwzzwDw6KOPcv755zN79myys7Pd+z/99NN8+umnTJo0iVmz8iko6OBuxcUtt9xCbW0t48aN47777iM/Px+AjIwMXnnlFS677DJyc3OZNWsWu3btIiUlhRtuuIGJEycyf/58pk+f7tP38emnn3Leeef5tG9QkVL6/MjPz5dm56GHtOml2Y8Ol06ns9V7nxR8IrkzW4KUL7wQpgGGGWCD9ON/7uujJ5wbNpuUERFSJiToU5C1x7vv18tXt7wqF/36Tfe2yWdtdb/OPWOn+3XeWdvdryfO3S4tEXYJUg6dtktGxTZKkHL0aK2fpUulvOkmbd8f/zjcn947/pwbmzdvPiil3ODrY/v27SH9LN7YuHGjvPLKKw1pq6FBym++kXLvXkOa8w+HQ8pvvpGN334rZ8yYIW02myHN7tixo902X8+PXmmpRMSfYHrOxHYzQ6dmT4V47XZUTxtV9B0aGsDhgDFjtL+Tk7XnG/5xJ1e+cyX/2vSVe98t21tcI1u3N7lfb97icL/etjkap0Ob53BwwxiaG6IhqoY9e7R+nnoK9IxWR8thChMwdepU5s6di8OAf0xYy7S4Oj1cVMSjjz6K1Rr+MHnvE5WjNhzxR5jaf2q795JjkhmZMZSohBp34FXRd9BjmgMH6s/aD7KmWvL7Mes4z/oUQkBMjCS+ahoAcXGSqBNatk5CgoSSSQCkpko4oRXwG+2RTDh0zjr368OH4YMPWvetMA/XXnute/Jjj8UlKqOGDOH0008P71hc9DpROXCkARKLyB+Q3+H7+dn5OOOKlaj0QfQL+6BB2nN8WhUAi4ZdxS8uO4n3VwjS0iAzU1BXp1m5o0YJmpu112PHtli+eXnaa6sVzjlH2zZ8OPzq8tMASBu3rdW8FWWp9F7CXlAy7J23pteJStFxIOE4+dmdi4o9tpDCIltoB6YIO/qFXReVyqjvABgWM829T1lZS/pvQgL0799yvKfbbMoU7fXcuS3bx42DMaO0WiwV6Su48bZK97HKUum9KDFpTa8SFSmhujyO+NRqshOzO9xnSvYUiC/haFFziEenCDf6hT09HdLSJIctn2Gx2tnxXVSr/XQLo7YWdO/IxIktc1XuvLNlkuQZZ7TsM326NhF8wuQG5Li3iJr6Bpdfrr2nLJXejxIXjfBHdQykogKcdisDB3TuJ52YORES3qK8sIf7UhV+o4uK1Qpvf7KfuX//X+LX38WGDdrP4N57NctjyhRcbjAYOlSbJPvii7B6tVbF+PbbYc8e+OIL+NGPNMulogLuuANiYuC7b2PI/P0RNpeu57XXbmbjRiUqvRkzBOrNJCq9ylI5VqRNJhgxOL7TfbLis4hNrqG+OoZmZaz0KfQLu9UKR1kHkY0MH+F0V2G46y6t0kJaGjz7rFY08oc/1Nayz86GSy/VrJfkZM0q2b1bW+0vJgZ+9SvtGbRZ2PnZ+Wwq2gRoloxyfxlPREQEeXl5TJw4kQsuuIDKykrvB3mhsrKStLQ0dwHHr7/+GiGEuzBlVVUVqampOD3q/ITrev7uu++yo6PZvG1YtmxZQKX8A6VXicrmvVqe8IThqZ3uI4RgYLbm7vCc7Kbo/egX9ogI2Fi0kRhrDJct1kpeZGZCm3JQ3SI/O5/tpdtptDcSEaEslWAQGxvL5s2b2bZtG6mpqTz33HPdbjMlJYXs7Gx2uma/rlmzhilTprBmjVaFeO3atZx00klYLO0vncESl87Snt9991126LN0u+j85ptv5qqrrgrG0DqkV4nKt3uPA5A/emCX+40YpJWzLi42j8moCD6e7q+NRRuZnDWZxRdrbtDRBtcYzR+Qj91pZ2vxVqxWZakEm1mzZnFUNzmB3/3ud0yfPp3c3Fz+53/+x739wQcfZMyYMZx88slcdtll/P73v2/X1uzZs90ismbNGu64445Wf+sFIF988UWmT5/O7NmT+eUvv099fT0Ab731FhMnTmTy5MmceuqpAGzfvp2TTjqJvLw8cnNz2bt3L6AVvtS333TTTW4BSUhI4Oc//zmTJ0/m66+/5p577nEvPvaLX/yCNWvW8N5773HXr35F3uWXs//IEfbv388555xDfn4+p5xyCrtctYaWLFni/pwdls43mF4VU9l18AQAJ48f0eV+E4an8x9gW0E5U6eq+hl9Bf2GLyJC8m3Rt1yZeyUjR8K558LMmcb2pWcfflv0LRERJ/VqS+Xaf107eFtJ16XvsUHcet9L3+f1z+Opc3wrfe9wOPj444+57rrrAG2Nkr1797J+/XqklCxYsIDPP/+c2NhY/vGPf7BlyxZsNhtTp051l0zxZM6cOaxevZrrr7+egoICFi9ezB/+8AdAE5V77rkH0KoF33DDDVRXwx13/IZ//vNl8vN/ygMPPMB//vMfBg4c6HbJLVu2jNtvv50rrriC5uZmHA4HO3fu5M033+Srr74iMjKSW2+9lddee42rrrqKuro6ZsyYwRNPPEF5eTnXXXcdu3btQghBZWUlKSkpLFiwgPPPPpuLR48Gq5V5N97IsmXLGDVqFOvWrePWW2/lk08+aff52pXOX7XK5/+LL/QqUTlY2IyIqmNQRnKX+00bMwDQ3GVXoUSlr6BbCw3OGmqaaxidppkn+gRFIxmcPJioiCgKThRgtSr3VzBoaGggLy+Po0ePMm7cOM466yxAE5WPPvqIKa6879raWvbu3UtNTQ0LFy4kJiaGmJgYdyn5tsyePZtHHnmEAwcOMHToUGJiYpBSUltby8aNG5kxYwagFYf8zW9+Q0VFJZWVtcyZMx/QROmaa67hBz/4gbtM/axZs3j44YcpLCzkoosuYtSoUXz88cds3LjRXduroaHBXQQzIiKC73//+wAkJycTExPDddddx/nnn8/555/fbsy1dXWsWbOGxYsXu7c1NTW12w86KJ1vML1KVI4fjSI2tQLoPFAPMGvsMAD2HKrpcj9F70IXlfIGLfaWk5wTtL4swsKQ5CEcqjrU6wP1yxd6L32/Iwil7/WYSn19PfPnz+e5557jtttuQ0rJr371K2666aZW+z/l46Jfo0aNorKykhUrVrjXN8nPz+dPf/oTQ4cOda8Gec011/Duu++SkzOZJ598hS1bPgM0q2TdunV88MEH5Ofns3HjRi6//HJmzJjBBx98wLnnnssf/vAHpJRcffXVPPLII+3GEBMT457tb7VaWb9+PR9//DFvv/02//d//9digbhiKU6nk5SUFJ/WhOmwdL6B9JqYipSSqqNZZA454XXfnLT+iLhyDh1REyD7Erq1UNbkEpWU4IkKaKKli4qyVIJHXFwczzzzDE888QR2u5358+ezfPly99K8R48epaSkhDlz5rBixQoaGxupra3l/fff77TNmTNn8vTTT7tFZdasWTz11FPueApATU0N2dnZ2Gw2PvzwNXesfP/+/cyYMYMHHniAjIwMjhw5QkFBAcOHD+e2225j4cKFbN26lXnz5vH2229T4poAVVFRwaFDh9qNpba2lqqqKs4991yefPJJtmzZAkBiYiI1rhLbSfHxDBs2jLfeegvQrof6fqGm11gqx6qLcZaNYPiZ3lPshBDE9Kug+Hiv0VSFD+g3ZWWNWkLH0JShQe0vJzmHlftWMt7asuCXIjhMmTKF3NxcXn/9dX74wx+yc+dOtyAkJCTw6quvMn36dBYsWEBubi5ZWVlMmjSJ5OSOXeVz5sxh5cqVTJumVVuYNWsWBQUFzJ49273Pgw8+yIwZM0hNzWDkyBk0NmoX+Lvuuou9e/cipWTevHlMnjyZxx57jL/+9a9ERkbSv39/7r33XlJTU3nooYc4++yzcTqdREZG8txzz5GT0/pmR3fbNTY2IqVk6dKlAFx66aXccN11PCMlbz/+OK+99hq33HILDz30EDabjUsvvZTJkycb/l17xZdSxvrDzOXN/7Z6rQQpf/bwNp/2z87bIqMGfxvkUZkP+nDp+08+0crQX/ToMzLhfxPaLY1gNPd/dr9kCXLemXY5c2ZQuzIEf86Nnlr6vqamRkopZV1dnczPz5cbN27sdpvl5Vrp+61bu92U/9TUaJ1v2GBos90pfd9rLJWvN2sLis/OS/Np/+xsJ0X706lrriM+qusYjKJ3oFsqJQ3HyMnMabc0gtHoMRsbjTgc6hwzAzfeeCM7duygsbGRq6++mqlT21czDxRVpkXDXKKyaxf89Kct9TL8YOsOzb9wSr5vi4HnDIpmU21/dpVuI39gnt9DVfQ89LhGccNRRgUzniIlLF9OzjFtbkPz4b3Ya/rDwy8H1t6QIdrUfkW30Zf1NRLTiImUWk2hMGMeUSkvh3nztHVZv/4a1q2DCRN8PvzAvigi4qrIyuw6nVhn7LAkkFa+2XtQiUofQbdUjtcd5czkMcHr6NgxuP56hqYAP4PG4oM4qwT85jeBt7lwISQlGTZEhXGYsPxWWDFPpPqNN7Qf4zvvaLeUrslGvlJyMIOUwUU+C3XuSG1+ytZ9qlZLX0EXlRp7BYOTBgevo0Zt1ciBv1sGQFPuKBzjc6G52f+HKyirIv3mJ6wFJcM2gPaYR1TeektbkGLRIm3Vo3/8AzyKtnVFs91G49GRDBlV5XN3w3O0XO2d++sCGq6i5+FO67XY6Z/Qv8t9u4VNS1WPjE8iLTaNZhqxOwRERvr/iItr1abCfJjK/WUCzCEqxcXw+eegzwa9+OIWN5gPrN1RCI39mDjJNxECTb8A9u2M9Xe0ih6Ke56XcJCVkBW8jnQBiIwkKyGLZmdd4PNUIiNbt6kwLSa5pocdc4jKp59q/xG9bML552sBJx9r0ny6rgyAmVMTfO4yMRGS+pdSXBDEO1aFqXCLisVOVnyIRCU+i0ZnfeAz6pWodIo+s92TJUuWMHDgQHdJ/Pfee6/V+1JK0tPTOXFCmyRdVFSEEIIvv/zSvU9GRgbl5eU+jyNYMZXPPvvMXciyq87fW72aR195xTSqZg5RWbdOW4xCn6iTnAzjx2vbfWDDZs3ffPasjld77IycsZXYjo6jvN73E0jRc/F0fwXVUvEoh5yVkEVjdywVq7V1mwqv3HHHHWzevJm33nqLa6+9ttXaJ0IIZs6cydcuL0jb0va7d+8mLS2NtDTfpiZA96/lnZW291VUFpx2Gvdcc033BmEg5hCV9eshP7/lrgxgxgxtuw//sd3boxHJhYwY4PuJAJA7yQkVI9lyZJ+/I1b0QFosFQeZ8b6lngdEW0vFUaPcX2Fg3LhxWK1WysrKWm3vqLS9p8jopVhWrFjBjBkzmDJlCmeeeSbFxVp5n9WrV5OXl0deXh5TpkyhtraGsrIirr/+VLeFpJeU/+ijj5g1axZTp05l8eLF7tIxQ4cO5e6772bq1Km89dZbPPPMM+7S9pdeeikHDx5k2bJlPPnkk+Tl5fHFF19QWlrK97//faZPn8706dP56quvAHhlxQp+8vjjICXXXHMNt912G7Nnz2b48OG8/fbbwf+i2xD+lGKbDTZtgltvbb19xgxYvhwKCmBE16Xsj+1LJynnAEIM8qvr2dMTeA0Ln64r44wgZpgqzIEuKsmxCURFRHW9c3fwFBVrFjbZiN0ugQDmEPQAUbn2WgZv24aXuvY57pwDX8jLAx/rP3bKunXrsFgsZGRktNo+Z84c7r//fgDWr1/P/fffz9NPPw1ooqKXYjn55JNZu3YtQgheeuklHn/8cZ544gl+//vf89xzzzFnzhxqa2uprIzhww//yMyZ83n++V/jcDior6+nrKyMhx56iFWrVhEfH89jjz3G0qVLue+++wBIS0tj0yZtddABAwZw4MABoqOj3aXtb775ZhISEvjFL34BwOWXX84dd9zBySefzOHDh5k/fz47XcICuG/Ai4qK+PLLL9m1axcLFizg4osv7t4X6SfhF5WtW7UUTFc5aTf63+vWdSkqjY1QVzSYCTP9L542/2TtbnX9N04I3cJoijChWwvpCQYu8dgRbSwVLPXY7E4gwv+2eoComI0nn3ySV199lcTERN588812lROmT5/Ot99+S11dHTabjYSEBIYPH86+fftYs2YNP//5zwEoLCzkkksuoaioiObmZoYN06qbz5kzhzvvvJMrrriCiy66CKt1EOPHT+fBB68lM9PGokWLyMvLY/Xq1ezYscNt+TQ3N7vrkQFccskl7te5ublcccUVLFq0iEWLFnX4uVatWtVq+eDq6mpqa9pXWl+0aBEWi4Xx48e7ratQEn5R0Us1ty2XMGECREVp719+eaeHb9zSCM4YJkzy/0c3fGgk1uTj7NzU+fLDit6DbqlkJgX5/91OVPZhswfoeO8BorJ8OT6Uvj9keOn7zrjjjjvcd/cdERcXx6hRo1i+fLm7TMvMmTNZuXIlJSUljBmjuS1++tOfcuedd7JgwQI+++wzlixZAsA999zDeeedx8qVK5kzZw5//et/mDr1VP74x885cuQDrrnmGu6880769evHWWedxeuvv97hOOLjW0r3fPDBB3z++eesWLGChx9+mO+++67d/k6nk7Vr1xITE9Oy0TOhwGWp6KXttU2hD96HP6aybRvExsLw4a23W61a3u+2bV0e/vEaV+bXNP9Tg4WAjLF7OL5rmN/HKnoeblFJ8C/25jdtUooRDpX9ZTJmz57NU0891aq0/dNPP83MmTPdlk1VVRUDB2pLk//5z392H7t//34mTZrE3XffzfTp09m3bxdFRYdITc3i+utv4Prrr2fTpk3MnDmTr776in37tJhtXV0de/bsaTcWp9PJkSNHmDt3Lo899hhVVVXU1ta2Km0PcPbZZ/Pss8+6/968ebOap9Ih27ZpVomlg6FMnOhVVNZubITIWk7JGxBQ96PzyrFVDODwEd/nuCh6Jrr7KysxhKISnwUWuwrUB4H6+noGDRrkfiy9916fv6c5c+ZQUFDgFpWpU6dSWFjYqrT9kiVLWLx4Mfn5+aSnaxU4HA5YsuQpJkyYSG5uLpGRkcyd+z02bvyMyy+fTH7+FN58801uv/127PYMnnzyFS677DJyc3OZNWuWe914TxwOB1deeSWTJk1iypQp3HbbbaSkpHDBBRfwzjvvuAP1zzzzDBs2bCA3N5fx48ezbNkyA75F4wm/+2vbNvje9zp+b+JEeO01qKrS0ow7YM+uCMjcxph03+uEeTJjloPVf4B/f3qCm64K8sVGEVYam+2Alf5JQV5C2iOlODM+EywOnI4AC/2plOJO8UwVprlZi88OGQKZ3jP7Fi9e3Mo1FB0d3W753YULF7Jw4cJW2+rr4fbbn2X4cEh1eVGPHoXzz7+a88+/milTwLVgIzt2wKRJZ/DNN9+0699zGd/IyMhW82R0Ro8ezdatW1tte/PNN1vvVFrKNRdcwDUXXABS8sorr7R6W882CyXhtVTKyuD4cU08OkLfvn17p02Ul0QRk1pGYnRiQEOYOUnLDNmypyKg4xU9h+qGeoDgi4qHpRJtjSYq0orTGaCoKEvFN3SBCbILSG/eU88669LpDIFHyiQuL0/CKyq6a8ubqHQQtNKprUggNSPwH1z+8OFgsbG/UK1X39upadTqvGUmhE5UAOKjo5FOS2C/fyUqvtHR1T4IeNOutiGOkF7zTSIw4RUVPT2us6yQIUMgPh527uzw7aYmsNclMyA78I8xKHkgIr6UwmPqR9vbqW5oAOEIg6ho2ToBxVWUqPhGiOrPd9RNV69DaqkoUUFblCshAVwZFu2wWGDMGG2/jg4/pLmshg/2Y1ZV2y6EhZiUKkpVBfxeT02jJioZ8Rned+4ObUQlIUZL8exlouJ0BuzTCwJhFJWu9g2y4RSUz9vdNOTwi8rYsV2vVjZuXKeWyvrdh7VdhvXr1jCS05qpKo/xvqOiR1PX2AgWOxlxIRaVaGKPoHwAABdYSURBVC3dPaBYu3lFZVtpaWmyaYTFJDEVz9chial4YkBnUkrKy8tbz4Xxk/Bmf+3aBaef3vU+Y8dqGWB1dZorzIOt+0sAyBvZvUrDWZmS4wf60WBrIDZSlcLvrdQ1NkKEnX6x3bsJ8Uo7S0WzpAMq1WJSUbHb7dcfP378pePHj0/Eh5vT8vLydjPbDaWxUUv8aWqCDmaZG0V9vdaNzQaVldq2ioqWLnfvbknYKynRMsEiAiik4DNVVS0D2bcPPCY+BkpMTAyDBvlX8sqT8IlKbS0cOaKJRlfo7+/ZA1OmtHprz0FtUa4pI/2rTtyWIQOj2VKbxZ7yPUzun+v9gIoK7azKCmKlW4Xh1DU1ISxOLCLIBrpHSjFAUqwmKpX1tSQn+5mlaNKU4vz8/BJgga/7T5s2TW7YsCF4A/rwQ21qwi23wPPPB62bv/0NrrgCfv1reOghbdvNN7csVFtQAK5qLuTmQv/+2mUuaDz8cMsy1atXw6mnBrEz3wif+2v3bu3Zm6joq2l14AI7dLQZgOz+3bsVGDk4CRwxfHvIS7Xiigq46CJIT9fOltNOg0OHutW3InTUNTVjiQjBJNc2lkqiy1IpqQ1giQWTWiqmQ59j0twc8m489V6Pmzmd2vYgD6d1oM4kNx7hExU9+O5NVEaO1OzHDkSl6LgkKrGaqG4WnJ0wVPOxb9p7tPOdamo0EXn/fbjnHnjkEa0u2axZ2uwnhelpaG4mIiIETm6bTUsycVWJSI7V3LalNQHMhVKi4hv61TvIV/GOuvG8ruuv9X9XSEUl4LINxhI+99euXZpYjBzZ9X7R0VpdMN2ycVFvq6e6PJaMtEYgqVtDGTRAU6VtB0o63+mGGzRh+/e/4ayztG3nnQezZ2vLH3/5ZZCdp4ru0tBsC82/yGZrtTZQcqy2QmFp7Qn/21Ki4hsmE5UQDceUohI+S2XnTq2kvS9mxtix7dKKd5XtgupB9M/u/p2nHhrZe6S64x0++ADefBPuv79FUAAmTdKcqWvXtjhVFaalsdlORChuozoVFWWpBA0lKkpU3OnEvjB2rBao9/jStpfsgPLRTBgb2cWBvtHflTx27KjA5mjz47Xb4Y47tDHcdVf7gy+7DM48U4vc6VkYCtPhcDpoarYRGQZRSYnTRKUsEEtFCM0CVqLSNSYVFbs9yHNVlKi4sNth717/RKWpCTyKsH2z5xA09mPG5I4LTfpDVhbEJTbjLBnN3oq9rd989VVtrI8+2rFVJQT87neaoHR3qTpF0DjReAKcViIjQzCtoo2o6DPqy+oCvOmIjFSi4g39Kt6mKGQouunouu75flB1TomKiwMHtG9az+zyhi4+Hi6wDd9p1TfHj+u+k1wIGD22GUomsr3Eo3ilw6Gl7E2ZAgu6yJ7My4MLL4Qnn4TqTlxoirBSXFsMTitR1hAEVdqIii5kAVkqWgNKVLxhUksl6EMyoaiEJ1Dva+aXjqeonHceAHt2Cb+a8EZ+XgybX53AluP/ZfGExdrGFSu0CUV//3vXs/4B7r0X3nkHli+Hn/3MmEH1JVatapeMYSTFzbtBziWqvhaee9X/BsaPh7lzfdvXbm+ZX0JL/kbAlorVapp0UdMSxpRiJSqtCY+o6IUkfVWE1FRtjQTXcZWNlZQfySAypplBg7qZT+xi8iQrNKRppV/muTYuXQo5OZoV4o1p0+Dkk+Hpp+GnP1WZYP5y/vlBdV0UTwScZxFTehx+8hP/G0hO9j1m1sZScYuKslSCh0ksFV37QyYqHU2SCTPhEZXt27Uikikpvh8zYYJ7XZUtx7dA2ViGDG/AYjFGVPQq+5u/c/2Ttm6FL77Q4iVWH7+mn/1MSy9euRIuuMCQcfUJbDZNUG67rWV2sMEUb/kD3GIlesJ4+G8XqeMd8dvfatl9Unq3WKGdqOinT7mKqQQPk4iKslTCKSoT/FypceJE+NOfQEo2HdsMR68g/0JjBEVvHqB052gqGipIfeEFiImBH/3I90YWLIABA7QyEUpUfKehQXvOyYGM4BR7LBH1CGklJi4GMvys75aToz3bbL6lwHdiqZyoq0FK6X8NLCUq3jFxoD6oQzKhqIQ+UO9waG4sf0VlwgStXtjhw3y2vhwa0jlvvnHFHzMy4KTTS+GbW1i3fZ2W9XXJJZDmxxLDkZFw441aHaKCAsPG1uup11ZkJC7wJQy8UVxbTKSIw2oNIPtLH5c+Tm90Yqk4HFLLQvMXJSreCbOlot8nhMVSadt5mAm9qBw4oFUU7Wy1x87wWFp4w1faDPozzjB2aL++JxLqM/nT44c0Abv5Zv8bue46rTzHSy8ZO7jeTChEpa6YSGIDC3V1U1TcfTqtWhaav/RlUWlo8C2Bwx9RKSmBoqKAhtNRN3Z7iwEbFlHROzdJMkfoRUVfbz4QSwWo3/Ytx7aOJXVQKd2oztwhF5ydgiWxhG82J2mz5WfM8L+RQYO0DLXly0MwnbaXECpREdE+h8daYZSoyAiK65So+MWLL2op/d58SP5kf910E1x1VUDD6Sz7yxSi0mctla1btefOlhDujJQUGDSI9bu/hIOnMn1OneFDEwIycw5y9MQEzY0V6PoPN90ExcVaSrLCOyFyf1mJCY2otEkpdr90RgRmqfTllOLjxzVrxdv8L38slePHtUcAdOb+UqLSQuhF5dtvYdQoSPRzXQmAvDz+uR9oTuL75/mROeYH42N3YqsYz+Fz5wfeyDnnaBbLH/9o3MB6M0EWFSklJXUlRMjowNxfsa7YnQHur5I6PzPPoG9bKrW1rZ87wx9Rqa313p4f3ShRaU14RGXq1MCOnTqVTyvyAFh0ThBEpaaG049+BM5I/rH5WODtRETA9dfDRx+pgL0vBFlUqpuqaXI0YSEqrO4vC1HK/eUv+pKK3lZz9EzL8rasbk1NwKtDdpz9JWkWNa7X7d8PevaXvtpjnxSVigqtfleAoiLz8thbdQbJA/cHJ/P01VdZVLIOgP9+Vda9tq67TruaqOrF3gmyqByr0W4QImR4REXvMzkqjaM1Aay905dFxV9LBby7Cg22VCrqqihv0s6xxmZbu/eVpRJMNm/WntssC+wre4an0VQ0i7yhe73v7C9SwgsvMGFiMhGxNWzc0M0Z8YMGafNWXn5Zy3ZTdE6QReVQlbY6p4WosGZ/ZcT053DVYf/7V6Lin6h4u4rX1mrmQwDfaVtRsTlsFFYVYY3SShF/uOcjv4fTLfq8pbJxo/YcoKi8sbsYbAmcG7nHwEG5+OIL+O47LD++haGTjlGycxQ1TYGZyG5uvRXKy7W1WBSdE2xRqdRERTgjw2qppMdmusfiF0pUfHd/tX3dFr16g2fbfuApKlLCuqPrsNslQ1K1RZnWH93k13C6TZ+3VNas0YL06ekBHb7yc80ldf6hD40clcZTT2kTHa+4gtNOjYDSCXywZW332pw3T8tye+op737evkwILBWrxYqQ1rDGVFKjszhSfQSn9HOBjb4sKrqYeBMAX2vNe7YTgKh4dmO3w3/2/QdkBGkJ2ty5/eUHqWysVKISEqTUltydMyegw53SyXeborFG1zL2wIdQWmrc2Pbvh3/9S0sjjo3lsnMHA/C3lQHcVXoiBNx+u+b2++yz7o+zt6JfrGONq5DgyeGqwwxKGoTdLsKSUqyLSlpMJs2OZo7X+pnO2pdTio12f3VTVNp281HBR0Rb4omL0f7f0mHhkwOfhLagZEvJhiB25DuhE5Xdu6GsTKvkGwDrCtfR8P/tnXlwFFUexz/dc2VyQQjBBF0DLGdEEIO4Qa1FEUSXKg7BWl0vtgRcFbZ0F0Qt1BVKVuRQUBZL8Sy0FAVlxfWoEhURFSIg4EUVSEgAuROSTGYy3W//eOlkApOh5wrHvE9VKqlJ93s93a/f9/1+7/d+r+wCunc9gI6AtWsTd21W0shJkwC4vMSD5qznizUCEa+FccstMsPyE08k4ELPUmpr5Wg8ZHSfSHZV7qKwTSHBYIzJo+MMKbbe+XYeaaFH7QJLZUslFvdXpHCr0HJiiAALrab88AHWV6zH68hoNBY8egYfbP+gdaO/HA75k3Ki8uWX8neMovLWxlWwrx9Dh3SQE1Nr1iTmuioqZKLKceOgoACQeSS7FB2icnsRW/dvja98r1duR/zRR7BhQwIu+CyktjapCx93Hd1FYdtCDMN+wulmOByyzcXp/mrrliGLVuCAbVJZVOy6vwKBpoebZEvFqmb51vcRCDx6k6hcmHcRy39cjq/OsHU5cRMqKqeJNdt6ovLxx7LT7t496lOFECz7cB8Ybq69NgMGDpTlJYKZM6Vr7v77m308dFAW7Clm2fcr46/jrrvknjDTp8df1tlIEkWl3qin4lhFo6USk6iAvL44LZW2bpmcVFkqNjGMpntuR1QyM5v+bokEiEpmpvRerNi2iqK8InThahSVfuf050jdEXYc3N14XKuJSkpZKoGAzNw7fHhMqU++qfiGik09cboMaegMHw5bt8rklPGwfbtM/Dh+PHTu3OxfQwZlgJHGi//bGL8LLDsbHnhA3gM1t3IiPl/SRKW8qhxTmPG5vyAuUbHqdGpp5KTl8OvRX6OrO1VFJfR+23F/2RGVON1fdX6DSlOuNVpftokbe9+IYWiNotItpye53lw+/uUzKo19oBlU1thsN7FgiYrTmWKi8vnn8gEOHx7T6S989wL6r0MoKRGy77H2Koknt5YQclMorzesBVFSIn9XbCvk812fx16Pxd13y305Jk06bczU04ba2qRN0m/+bTMARXlFsbu/ICGiYhjyOr7f/310daeqqIR2+nYsFSv1U5IslW/L12MEHaRlyGfxz0sfYsrAKc2WimjCySe3fEKnrO5keJ3gCPD2lgR4O1oiZS2Vt9+WncbVV0d96sHagyx9fxfm3r6MHtXQI3TrBr16wbJlsV/TW29Jy+HRRyE//4R/5+dDryITfcttLPzmmdjrsfB6Yf58aWHNnRt/eWcTSXR/le4pxaE56Jvft3XcX4YhByxh91OB/h37s2nfJoJmFAOLVBWVaATA70+6qMz8bDYAReedD8BfLhiHx+k5Iaq3X0E/is8ZSMe2eaSlaew6vJd1u9dFVZdtUlJUjh2D11+HG26IqeOYs3YudR8/RG5ekAkTQv4xbpyc/LdS6UdDebmc5xgwQForLXD/VB1zbx+Wv+ePf8IeYORIuP56ePjhpuwCiuSKyt5SivKKSHelt477y7JCw4QUB4NQXFBMbX0tPx38yX7dqRpSHNrpR+P+SkL0108Hf+K/P8j1cW2yHM2qCbdUxO+Xn2V63XhowxNrkxT9mZKisnRpzBtelVWWMW+OC3YNYuZjzub9zu23y6e2aFF0hfp8MGqUbISvvhpx6HrTTdC5i4G26jkmvzk7/rkVTYPFi+U2kyNHJnatzZlMkkRFCEHp3lKKOxYjBJhmK1gqlkXRgqVS3LEYkBaUbVLdUklPj2xVmKYU3Wgm6r3eqCyVZ759Bhey/OOraSlLsdsNHrdOz5w+rPx5ZfRzaXZIOVGpqZHRVQMGRL3hlSlMxsydT/0njzB8dDUTJx53QF6e3Gjn+eflhLsdfD7ZmZeWyu2Ce/SIeLjLBe+ucOAOtmf1v+/hpdKlUX2HsLRvDytWyP1Whg6VaVxSnSSJyu6q3eyv2U9xQXHj+3YqRCXUUumR24MMVwbr96y3X7fLJV1qp0mn0WpYlkR+fmSrwrrndifq09Nl8IxNS6XKX8Urm19hRNexYasxjKbHfbyouN3w++xe6JrOovVRDoDtkHKi8sgjch3IvHlRRX0JIZi8Yjrrn7mH3IJqXn85M/zpM2bIRSUTJ558JFdWBoMGyVDkJUtgxAhb19KnD7z2sgP2DGDC3w/z1e6vbH+PFrnkEnj3XfjhBxkRsDUBrrUzmSSJyhtb3gBgcOfBjd6jmN1fXm/comIY4NAdDO4ymHd+fId6w6b1YZWVataKZUkUFES2Kqze3e6cSlaWVAablspTXz9FdaCaW3vfEbYaKwBE05qLiscjRcVFBqN7jWbxhsWx7acTiZQSlUWL5IT0nXdGlZrliO8IN70xkWenlaBXdeLdN7Nb3s8rPx8WLIDVq+Hmm8OPPI4cgccfl9sRb9sGy5fL+ZgoGDvWwR131mJ8NZlBf3uHlza+HL8r7Jpr4NNPobISioth2rTUdYclQVRMYfJc6XP8sfCP9MrrFW6qIzrisFT0hrfMeucnXDyBfdX7eO/n9+zVrUTFnqjYdX9lZkplsCEq5VXlzF47mzFFY+iZ0ydsNeH69VBLJRCAGVfOoLa+lumfJnit2mkoKrG+YuHx+WT6lIULYeVKuO462elHQAjBId8h1m7fwosrdvDhKjeBbY+Drz0LnjW5/PKTWDi33y4742nTpLiMGgXnny+3H920SYYz+/0ynHnBghPWo9jlPwvT2bHDz6er5vLX70p56A9PM+aajtxybQ8u7NiDNGda9IVedpm0Uu67D2bPlhbdsGEwZIgUmq5dpZsv1m2NzxQSLCqGaTB99XR2Ht3JrMGz5Gen0P2lac0XPA/rOoxObTsx9ZOpXFxwMV1yukQuM1VFJdT9VV0t5070MONga8bcrvsrM1P+nMT99cuhXxj95mgAZg2eRd2eE6sRQv6cTFR6tO/B5EsnM//r+ZyXfR5TL5uKx+mxcxcicxqKihbNiLt///5iQ0iqkfwuH+LzZzXdWassTWuy/TRAgECTv4WGaWqYhoYRcBH0exC+bKhrB4DTW8NVQwPMeDCHAQOi+Cbr1skcXqtXw9Gj8kXs2VNmCr7tNrjooigKC49pwpIlJjOePMru7e0avmsQMn9Dd/txevw4XEF0h0DTQNNE4+2Qf1j3OoxIGIZsfcGgrCgUWVjT3zZZ+OBRbr37T8cVpZUKIfrbLsQmx7eNf01dyrylneydbAX6p8UgzI0IENJCqTPqCBpB8jML6JbbFdA4fFgaqk8/HTHgr2WmTIE5c+CKKyIft3OnjC587TVpPTeg6/L1sE4/Fqhiy29bMUQQrzMdp+5AQ2toGsc944AffHWykAQNMM5pt4dftoxt9lmy2gac2D66FK3gUGWHyCdZ0RUejxSOlnyX1nFpaU17F7V0rCVMlq8qVKQabq0J1GsCv26iC+hdk0GboJNttZ04HGzDk50XMWXnXQCUZG1l3bHePFb4AjPLbiUg3FyRvZk1VX0ZlfsFewPt+PpYb67I3oxA8HOGjwOuenSh4TV1dEKedizODytKwOo3YvbvNufhCRX845Ebmn1mt33EZanoDoHuMK0aQW9Y2Wk5GEPukoYADXRNoOlWOqU6MjKraNf2KN07HWLM0HO5elAGbndG9BdTUiJdW0LIztnlCj+qiQNdh/HjdcaPb8eBA/DOh/v5aM1Byir8VNUEqfNBwK9jGLrUVyEa74AQoZ1BmNbj0MDtATyycRhG05qHZqItrAJtXO+ps3Acuoau2xw56YBHB0eUKeHD4NIcZDlyyE3PJS89D+uV7dABzj0XrroqxoJHjJBbYR8v+MfTrRtceKFMJRTCvffK0y1ynNlc8rt+7KveR02gmqBpENJampfp1sHQINqU+RFw6IkrKxZ03bTXPtwO8DjkLRERjndqsg0JveEZtXCs3lCmpoE/eMJxGuBEI11oZAZcFPjduIUOuknfzB2kO+oYm7+GDTU92R/IAeDqnA0Ma1/KMTOTDcdk8M+VbTfy5/zVHKnPwutosjB7+zwcDjo45KzHp5uYWtNTDzfWPCnWu+N0QH2E7x0lWhyDl7gsFcWZR2tZKoozj9a0VBRnHnbbR+tu0qVQKBSKsxolKgqFQqFIGFG5vzRNOwDEuR2i4hRTKITIS3Shqm2cFSSlbYBqH2cJttpHVKKiUCgUCkUklPtLoVAoFAlDiYpCoVAoEoYSFYVCoVAkDCUqCoVCoUgYSlQUCoVCkTCUqCgUCoUiYShRUSgUCkXCUKKiUCgUioShREWhUCgUCeP/GEUvkGOW7AcAAAAASUVORK5CYII=\n",
+      "text/plain": [
+       "<Figure size 432x288 with 6 Axes>"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "# Author: Remi Flamary <remi.flamary@unice.fr>\n",
+    "#\n",
+    "# License: MIT License\n",
+    "\n",
+    "import numpy as np\n",
+    "import matplotlib.pylab as pl\n",
+    "import ot\n",
+    "# necessary for 3d plot even if not used\n",
+    "from mpl_toolkits.mplot3d import Axes3D  # noqa\n",
+    "from matplotlib.collections import PolyCollection  # noqa\n",
+    "\n",
+    "#import ot.lp.cvx as cvx\n",
+    "\n",
+    "#\n",
+    "# Generate data\n",
+    "# -------------\n",
+    "\n",
+    "#%% parameters\n",
+    "\n",
+    "problems = []\n",
+    "\n",
+    "n = 100  # nb bins\n",
+    "\n",
+    "# bin positions\n",
+    "x = np.arange(n, dtype=np.float64)\n",
+    "\n",
+    "# Gaussian distributions\n",
+    "# Gaussian distributions\n",
+    "a1 = ot.datasets.make_1D_gauss(n, m=20, s=5)  # m= mean, s= std\n",
+    "a2 = ot.datasets.make_1D_gauss(n, m=60, s=8)\n",
+    "\n",
+    "# creating matrix A containing all distributions\n",
+    "A = np.vstack((a1, a2)).T\n",
+    "n_distributions = A.shape[1]\n",
+    "\n",
+    "# loss matrix + normalization\n",
+    "M = ot.utils.dist0(n)\n",
+    "M /= M.max()\n",
+    "\n",
+    "#\n",
+    "# Plot data\n",
+    "# ---------\n",
+    "\n",
+    "#%% plot the distributions\n",
+    "\n",
+    "pl.figure(1, figsize=(6.4, 3))\n",
+    "for i in range(n_distributions):\n",
+    "    pl.plot(x, A[:, i])\n",
+    "pl.title('Distributions')\n",
+    "pl.tight_layout()\n",
+    "\n",
+    "#\n",
+    "# Barycenter computation\n",
+    "# ----------------------\n",
+    "\n",
+    "#%% barycenter computation\n",
+    "\n",
+    "alpha = 0.5  # 0<=alpha<=1\n",
+    "weights = np.array([1 - alpha, alpha])\n",
+    "\n",
+    "# l2bary\n",
+    "bary_l2 = A.dot(weights)\n",
+    "\n",
+    "# wasserstein\n",
+    "reg = 1e-3\n",
+    "ot.tic()\n",
+    "bary_wass = ot.bregman.barycenter(A, M, reg, weights)\n",
+    "ot.toc()\n",
+    "\n",
+    "\n",
+    "ot.tic()\n",
+    "bary_wass2 = ot.lp.barycenter(A, M, weights, solver='interior-point', verbose=True)\n",
+    "ot.toc()\n",
+    "\n",
+    "pl.figure(2)\n",
+    "pl.clf()\n",
+    "pl.subplot(2, 1, 1)\n",
+    "for i in range(n_distributions):\n",
+    "    pl.plot(x, A[:, i])\n",
+    "pl.title('Distributions')\n",
+    "\n",
+    "pl.subplot(2, 1, 2)\n",
+    "pl.plot(x, bary_l2, 'r', label='l2')\n",
+    "pl.plot(x, bary_wass, 'g', label='Reg Wasserstein')\n",
+    "pl.plot(x, bary_wass2, 'b', label='LP Wasserstein')\n",
+    "pl.legend()\n",
+    "pl.title('Barycenters')\n",
+    "pl.tight_layout()\n",
+    "\n",
+    "problems.append([A, [bary_l2, bary_wass, bary_wass2]])\n",
+    "\n",
+    "#%% parameters\n",
+    "\n",
+    "a1 = 1.0 * (x > 10) * (x < 50)\n",
+    "a2 = 1.0 * (x > 60) * (x < 80)\n",
+    "\n",
+    "a1 /= a1.sum()\n",
+    "a2 /= a2.sum()\n",
+    "\n",
+    "# creating matrix A containing all distributions\n",
+    "A = np.vstack((a1, a2)).T\n",
+    "n_distributions = A.shape[1]\n",
+    "\n",
+    "# loss matrix + normalization\n",
+    "M = ot.utils.dist0(n)\n",
+    "M /= M.max()\n",
+    "\n",
+    "\n",
+    "#%% plot the distributions\n",
+    "\n",
+    "pl.figure(1, figsize=(6.4, 3))\n",
+    "for i in range(n_distributions):\n",
+    "    pl.plot(x, A[:, i])\n",
+    "pl.title('Distributions')\n",
+    "pl.tight_layout()\n",
+    "\n",
+    "#\n",
+    "# Barycenter computation\n",
+    "# ----------------------\n",
+    "\n",
+    "#%% barycenter computation\n",
+    "\n",
+    "alpha = 0.5  # 0<=alpha<=1\n",
+    "weights = np.array([1 - alpha, alpha])\n",
+    "\n",
+    "# l2bary\n",
+    "bary_l2 = A.dot(weights)\n",
+    "\n",
+    "# wasserstein\n",
+    "reg = 1e-3\n",
+    "ot.tic()\n",
+    "bary_wass = ot.bregman.barycenter(A, M, reg, weights)\n",
+    "ot.toc()\n",
+    "\n",
+    "\n",
+    "ot.tic()\n",
+    "bary_wass2 = ot.lp.barycenter(A, M, weights, solver='interior-point', verbose=True)\n",
+    "ot.toc()\n",
+    "\n",
+    "\n",
+    "problems.append([A, [bary_l2, bary_wass, bary_wass2]])\n",
+    "\n",
+    "pl.figure(2)\n",
+    "pl.clf()\n",
+    "pl.subplot(2, 1, 1)\n",
+    "for i in range(n_distributions):\n",
+    "    pl.plot(x, A[:, i])\n",
+    "pl.title('Distributions')\n",
+    "\n",
+    "pl.subplot(2, 1, 2)\n",
+    "pl.plot(x, bary_l2, 'r', label='l2')\n",
+    "pl.plot(x, bary_wass, 'g', label='Reg Wasserstein')\n",
+    "pl.plot(x, bary_wass2, 'b', label='LP Wasserstein')\n",
+    "pl.legend()\n",
+    "pl.title('Barycenters')\n",
+    "pl.tight_layout()\n",
+    "\n",
+    "#%% parameters\n",
+    "\n",
+    "a1 = np.zeros(n)\n",
+    "a2 = np.zeros(n)\n",
+    "\n",
+    "a1[10] = .25\n",
+    "a1[20] = .5\n",
+    "a1[30] = .25\n",
+    "a2[80] = 1\n",
+    "\n",
+    "\n",
+    "a1 /= a1.sum()\n",
+    "a2 /= a2.sum()\n",
+    "\n",
+    "# creating matrix A containing all distributions\n",
+    "A = np.vstack((a1, a2)).T\n",
+    "n_distributions = A.shape[1]\n",
+    "\n",
+    "# loss matrix + normalization\n",
+    "M = ot.utils.dist0(n)\n",
+    "M /= M.max()\n",
+    "\n",
+    "\n",
+    "#%% plot the distributions\n",
+    "\n",
+    "pl.figure(1, figsize=(6.4, 3))\n",
+    "for i in range(n_distributions):\n",
+    "    pl.plot(x, A[:, i])\n",
+    "pl.title('Distributions')\n",
+    "pl.tight_layout()\n",
+    "\n",
+    "#\n",
+    "# Barycenter computation\n",
+    "# ----------------------\n",
+    "\n",
+    "#%% barycenter computation\n",
+    "\n",
+    "alpha = 0.5  # 0<=alpha<=1\n",
+    "weights = np.array([1 - alpha, alpha])\n",
+    "\n",
+    "# l2bary\n",
+    "bary_l2 = A.dot(weights)\n",
+    "\n",
+    "# wasserstein\n",
+    "reg = 1e-3\n",
+    "ot.tic()\n",
+    "bary_wass = ot.bregman.barycenter(A, M, reg, weights)\n",
+    "ot.toc()\n",
+    "\n",
+    "\n",
+    "ot.tic()\n",
+    "bary_wass2 = ot.lp.barycenter(A, M, weights, solver='interior-point', verbose=True)\n",
+    "ot.toc()\n",
+    "\n",
+    "\n",
+    "problems.append([A, [bary_l2, bary_wass, bary_wass2]])\n",
+    "\n",
+    "pl.figure(2)\n",
+    "pl.clf()\n",
+    "pl.subplot(2, 1, 1)\n",
+    "for i in range(n_distributions):\n",
+    "    pl.plot(x, A[:, i])\n",
+    "pl.title('Distributions')\n",
+    "\n",
+    "pl.subplot(2, 1, 2)\n",
+    "pl.plot(x, bary_l2, 'r', label='l2')\n",
+    "pl.plot(x, bary_wass, 'g', label='Reg Wasserstein')\n",
+    "pl.plot(x, bary_wass2, 'b', label='LP Wasserstein')\n",
+    "pl.legend()\n",
+    "pl.title('Barycenters')\n",
+    "pl.tight_layout()\n",
+    "\n",
+    "\n",
+    "#\n",
+    "# Final figure\n",
+    "# ------------\n",
+    "#\n",
+    "\n",
+    "#%% plot\n",
+    "\n",
+    "nbm = len(problems)\n",
+    "nbm2 = (nbm // 2)\n",
+    "\n",
+    "\n",
+    "pl.figure(2, (20, 6))\n",
+    "pl.clf()\n",
+    "\n",
+    "for i in range(nbm):\n",
+    "\n",
+    "    A = problems[i][0]\n",
+    "    bary_l2 = problems[i][1][0]\n",
+    "    bary_wass = problems[i][1][1]\n",
+    "    bary_wass2 = problems[i][1][2]\n",
+    "\n",
+    "    pl.subplot(2, nbm, 1 + i)\n",
+    "    for j in range(n_distributions):\n",
+    "        pl.plot(x, A[:, j])\n",
+    "    if i == nbm2:\n",
+    "        pl.title('Distributions')\n",
+    "    pl.xticks(())\n",
+    "    pl.yticks(())\n",
+    "\n",
+    "    pl.subplot(2, nbm, 1 + i + nbm)\n",
+    "\n",
+    "    pl.plot(x, bary_l2, 'r', label='L2 (Euclidean)')\n",
+    "    pl.plot(x, bary_wass, 'g', label='Reg Wasserstein')\n",
+    "    pl.plot(x, bary_wass2, 'b', label='LP Wasserstein')\n",
+    "    if i == nbm - 1:\n",
+    "        pl.legend()\n",
+    "    if i == nbm2:\n",
+    "        pl.title('Barycenters')\n",
+    "\n",
+    "    pl.xticks(())\n",
+    "    pl.yticks(())"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.6.5"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 0
+}