summaryrefslogtreecommitdiff
path: root/notebooks/plot_OT_L1_vs_L2.ipynb
blob: d6a87619965888ce9cbe3748d490cde436a75943 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "# 2D Optimal transport for different metrics\n",
    "\n",
    "\n",
    "2D OT on empirical distributio  with different gound metric.\n",
    "\n",
    "Stole the figure idea from Fig. 1 and 2 in\n",
    "https://arxiv.org/pdf/1706.07650.pdf\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "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",
    "import ot.plot"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Dataset 1 : uniform sampling\n",
    "----------------------------\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAADSCAYAAAAffFTTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAVU0lEQVR4nO3df5RtZX3f8fcngL9AC8gtooJXhZjQVNE1IFmLNhA1AeMqJjUUWlz4qyRUW21pbok/4IqxWbltiU3baLEiBAx6o0Zpgl0QHIOmBpiroCghIgGBXC4jBAG1yI9v/9j7yuHcmTvD7L1nztx5v9Y665yz9z77ec4zZ873PHt/n/2kqpAkadL8xEpXQJKkuRigJEkTyQAlSZpIBihJ0kQyQEmSJpIBSpI0kQxQWvOSvCHJl1a6Hn1KUkkObh9/KMl7etrvQUkeSLJb+/wLSd7Sx77b/X0uySl97U+rmwFKJDkqyf9N8r0k9yT5iySHr3S9JkGS9e2X/e4rWIdbkrxyqa+vql+vqvf1UU5Vfaeq9qqqR5Zan5HyNia5aGz/x1XVBV33rV3Div3TaTIkeQbwJ8BpwGbgScA/Ah4coKzdq+rhvvc7yXal97wrvRetDvag9JMAVXVxVT1SVT+sqsuq6msASX4iybuT3JrkriR/kOTvteuOTnL76M5Gf4W3v5A/meSiJPcBb0iyW5J3Jvl2kvuTbElyYLv9TyW5vO3F3ZjkhPkqneSNSW5o93Fzkl8bWXd0ktuTnN7WeWuSN46sf2aSS5Lcl+Rq4IU7aZ8r2/t720NbP5vkhUk+n+TuJN9N8rEke4+1wX9I8jXg+0l2T/KyJF9t6/tHST6R5LdGXvOaJNcmubftzb64XX4hcBDwv9vyN8zTHr/Rvs+/TfKmsXXnby8ryX5J/qQt554kX2z/xjuUM9J7fHOS7wCfn6dH+cIkV7ft+dkk+47+HcbqckuSVyY5Fngn8M/a8q5r1//4kOECn73t9TglyXfav8O7Rso5IslMW6dtSc7Zyd9Yk6qqvK3hG/AM4G7gAuA4YJ+x9W8CbgJeAOwFfBq4sF13NHD72Pa3AK9sH28EHgJeS/Nj6KnAbwBfB14EBHgJ8ExgT+A24I00PfuXAt8FDp2n3r9EE1gC/BzwA+BlI/V6GDgb2AN4dbt+n3b9x2l6i3sCPwPcAXxpnnLWAwXsPrLsYOBVwJOBdTRB7ANjbXAtcGD7np8E3Aq8va3PrwA/An6r3f6lwF3Ay4HdgFPafTx5vE3nqeOxwLb2vewJ/GFb54Pb9eePlPXbwIfaeuxB01vOXOWMvPc/aPf71PH2AL7Qtt/2sj8FXPQEPh8Xja3/AvCWRXz2ttfjw229XkLT6//pdv2Xgde3j/cCjlzp/zVvT/xmD2qNq6r7gKN47J99tu1d7N9u8i+Ac6rq5qp6APhN4MQs/pzMl6vqM1X1aFX9EHgL8O6qurEa11XV3cBrgFuq6qNV9XBVfZXmy+5X56n3n1bVt9t9/DlwGc2X7XYPAWdX1UNVdSnwAPCiNCf3/ylwZlV9v6qupwnOi1ZVN1XV5VX1YFXNAufQBMlRv1dVt7Xv+UiaoPt7bX0+DVw9su2pwP+sqquq6cVeQPNle+Qiq3QC8NGqur6qvk/zxT+fh4ADgOe1dfliVS10Qc6NbVv9cJ71F46U/R7ghLadu1rMZ++91fT6rwOuowlU0LzPg5PsV1UPVNVf9lAfLTMDlKiqG6rqDVX1XJpfws8GPtCufjbNr//tbqX5st2fxblt7PmBwLfn2O55wMvbQ0/3JrmX5gvqWXPtNMlxSf6yPUx1L00vab+RTe6ux58v+QHNL+l1bf1H6zX6/haUZP8kH09yR3vo8qKxshnb/7OBO8YCwej65wGnj733A9vXLcazWfz7+U80vZLL2kOjZyxi/+N/w52tv5WmZzbeHkuxmM/enSOPt/+NAd5Mc/j6r5Jck+Q1PdRHy8wApcepqr+iOST0M+2iv6X5At3uIJrDZ9uA7wNP276i/dW8bnyXY89vY+5zPrcBf15Ve4/c9qqq08Y3TPJkmt7Vfwb2r6q9gUtpDvctZLat/4Fj72k+c/Uu/mO7/B9W1TOAk+coe/R1W4HnJBndZrT824D3j733p1XVxTupw6itLPL9VNX9VXV6Vb0A+CfAv0vyigXKWaj88bIfojk8u9DnY6H97uyzt1NV9a2qOgn4+8DvAJ9MsudCr9NkMUCtcWkSE05P8tz2+YHAScD2QyIXA/82yfOT7EXz5fyJtnfy18BTkvxSkj2Ad9Ocl9mZ/wW8L8khabw4yTNpMgl/Msnrk+zR3g5P8tNz7ONJbTmzwMNJjgN+YTHvt5r06E8DG5M8LcmhNOd85jMLPEpzHmS7p9McMvxekufQnFfbmS8DjwBvaxMmjgeOGFn/YeDXk7y8bZM92zZ9ert+21j54zbTJKAcmuRpwFnzbZgmGePgNlh+r63Xo4ssZz4nj5R9NvDJtp0X+nxsA9Ynme97aGefvZ1KcnKSdVX1KHBvu/jRnb1Gk8cApftpTs5fleT7NIHpeuD0dv15wIU0iQB/A/w/4F8DVNX3gH9FE3TuoPnF/LisrTmcQ/OFehlwH/AR4KlVdT9NkDmR5pfznTS/fHcIeO22/6bdz98B/xy45Am857fRHAq6k6a3+NH5NqyqHwDvB/6iPfx2JPBe4GU0X/B/ShPw5lVVP6JJjHgzzZflyTQB+cF2/QzwL4H/3r6fm4A3jOzit4F3t+X/+zn2/zmaQ7Kfb1/7+Z1U5xDgz2gC7JeB36+q6cWUsxMX0rTjncBTaP42i/l8/FF7f3eSr8yx33k/e4twLPCNJA8A/xU4cSfn0DShtmfvSFpGSa4CPlRV8wZHaa2zByUtgyQ/l+RZ7SG+U4AXA/9npeslTTKvJCEtjxfx2Nirm4HXVdXWla2SNNk8xCdJmkge4pMkTSQDlCRpIi3rOaj99tuv1q9fv5xFSpIm3JYtW75bVeOD/Jc3QK1fv56ZmZnlLFKSNOGSzHl5Lg/xSZImkgFKkjSRFgxQSQ5MMp3km0m+keTt7fJ900wu9632fp/hq6uJtGkTTE8/ftn0dLNckpZoMT2oh4HTq+pQmvlp3tpeYPMM4IqqOgS4on2utejww+GEEx4LUtPTzfPDD1/Zekla1RYMUFW1taq+0j6+H7gBeA5wPI9N9HYBzaypWouOOQY2b26C0plnNvebNzfLJWmJntA5qCTraaanvopmHp7tl2q5k3kmsEtyapKZJDOzs7MdqqqJdswxcNpp8L73NfcGJ0kdLTpAtfOxfAp4RztN+I+1M4XOec2kqjq3qqaqamrduh3S3LWrmJ6GD34Q3vOe5n78nJQkPUGLClDtZGOfAj5WVdvnvtmW5IB2/QHAXcNUURNv+zmnzZvh7LMfO9xnkJLUwWKy+EIzqdwNVXXOyKpLeGwm0lOAz/ZfPa0K11zz+HNO289JXXPNytZL0qq24NXMkxwFfBH4Oo9NmfxOmvNQm4GDgFuBE6rqnp3ta2pqqryShCRpVJItVTU1vnzBSx1V1ZeAzLP6FV0rpl3Apk1NSvloYsT0dNOD2rBh5eolaVXzShLqznFQkgbgjLrqbnQc1GmnNVl8joOS1JE9KPXDcVCSemaAUj8cByWpZwYodec4KEkDMECpO8dBSRqAAUqSNJEMUOrONHNJAzDNXN2ZZi5pAPag1A/TzCX1zAClfphmLqlnBih1Z5q5pAEYoNSdaeaSBmCAkiRNJAOUujPNXNIATDNXd6aZSxqAPSj1wzRzST0zQKkfpplL6pkBSt2ZZi5pAAYodWeauaQBpKqWrbCpqamamZlZtvIkSZMvyZaqmhpfbg9K3W3atOPhvOnpZrkkLZEBSt05DkrSABwHpe4cByVpAPag1A/HQUnqmQFK/XAclKSeGaDUneOgJA3AAKXuHAclaQALBqgk5yW5K8n1I8s2JrkjybXt7dXDVlOStNYspgd1PnDsHMt/t6oOa2+X9lstrSqmmUsawIIBqqquBO5ZhrpotRpNMz/zzMfOR5nJJ6mDLueg3pbka+0hwH3m2yjJqUlmkszMzs52KE4TzTRzST1baoD6IPBC4DBgK/Bf5tuwqs6tqqmqmlq3bt0Si9PEM81cUs+WFKCqaltVPVJVjwIfBo7ot1paVUwzlzSAJQWoJAeMPP1l4Pr5ttUaYJq5pAEsON1GkouBo4H9gG3AWe3zw4ACbgF+raq2LlSY021IksbNN93GgheLraqT5lj8kV5qpV3Dpk1NSvloYsT0dNOD2rBh5eolaVXzShLqznFQkgbgdBvqzuk2JA3AHpT64TgoST0zQKkfjoOS1DMDlLpzHJSkARig1J3joCQNwAAlSZpIBih1Z5q5pAGYZq7uTDOXNAB7UOqHaeaSemaAUj9MM5fUMwOUujPNXNIADFDqzjRzSQMwQEmSJpIBSt2ZZi5pAKaZqzvTzCUNwB6U+mGauaSeGaDUD9PMJfXMAKXuTDOXNAADlLozzVzSAFJVy1bY1NRUzczMLFt5kqTJl2RLVU2NL7cHpe42bdrxcN70dLNckpbIAKXuHAclaQCOg1J3joOSNAB7UOqH46Ak9cwApX44DkpSzwxQ6s5xUJIGYIBSd46DkjSABQNUkvOS3JXk+pFl+ya5PMm32vt9hq2mJGmtWUwP6nzg2LFlZwBXVNUhwBXtc61VpplLGsCCAaqqrgTuGVt8PHBB+/gC4LU910uryWia+ZlnPnY+ykw+SR0s9RzU/lW1tX18J7D/fBsmOTXJTJKZ2dnZJRaniWeauaSedU6SqOZifvNe0K+qzq2qqaqaWrduXdfiNKlMM5fUs6UGqG1JDgBo7+/qr0padUwzlzSApQaoS4BT2senAJ/tpzpalUwzlzSABafbSHIxcDSwH7ANOAv4DLAZOAi4FTihqsYTKXbgdBuSpHHzTbex4MViq+qkeVa9onOttGvYtKlJKR9NjJiebnpQGzasXL0krWpeSULdOQ5K0gCcbkPdOd2GpAHYg1I/HAclqWcGKPXDcVCSemaAUneOg5I0AAOUunMclKQBGKAkSRPJAKXuTDOXNADTzNWdaeaSBmAPSv0wzVxSzwxQ6odp5pJ6ZoBSd6aZSxqAAUrdmWYuaQAGKEnSRDJAqTvTzCUNwDRzdWeauaQB2INSP0wzl9QzA5T6YZq5pJ4ZoNSdaeaSBmCAUnemmUsaQKpq2QqbmpqqmZmZZStPkjT5kmypqqnx5fag1N2mTTsezpuebpZL0hIZoNSd46AkDcBxUOrOcVCSBmAPSv1wHJSknhmg1A/HQUnqmQFK3TkOStIADFDqznFQkgbQKUkiyS3A/cAjwMNz5bFLkrQUffSgjqmqwwxOa5hp5pIGYJq5ujPNXNIAuvagCrgsyZYkp861QZJTk8wkmZmdne1YnCaWaeaSetY1QB1VVS8DjgPemuQfj29QVedW1VRVTa1bt65jcZpYpplL6lmnAFVVd7T3dwF/DBzRR6W0yphmLmkASw5QSfZM8vTtj4FfAK7vq2JaRUwzlzSAJU+3keQFNL0maJIt/rCq3r+z1zjdhiRp3HzTbSw5i6+qbgZe0qlW2jVs2tSklI8mRkxPNz2oDRtWrl6SVjWvJKHuHAclaQCOg1J3joOSNAB7UOqH46Ak9cwApX44DkpSzwxQ6s5xUJIGYIBSd46DkjQAA5QkaSIZoNSdaeaSBmCaubozzVzSAOxBqR+mmUvqmQFK/TDNXFLPDFDqzjRzSQMwQKk708wlDcAAJUmaSAYodWeauaQBmGau7kwzlzQAe1Dqh2nmknpmgFI/TDOX1DMDlLozzVzSAAxQ6s40c0kDSFUtW2FTU1M1MzOzbOVJkiZfki1VNTW+3B6Uutu0acfDedPTzXJJWiIDlLpzHJSkATgOSt05DkrSAOxBqR+Og5LUMwOU+uE4KEk9M0CpO8dBSRqAAUrdOQ5K0gA6Bagkxya5MclNSc7oq1JaZTZs2PGc0zHHNMvHbNw4/26Wum6lXmudVna/2vUteaBukt2AvwZeBdwOXAOcVFXfnO81DtRVAvN95Ja6bqVea51Wvk7aNQwxUPcI4KaqurmqfgR8HDi+w/4kSfqxLgHqOcBtI89vb5c9TpJTk8wkmZmdne1QnFarjRubX8JJ83z7440bl76uy36t0+quk9aOLof4XgccW1VvaZ+/Hnh5Vb1tvtd4iE+TeJjIOq3eOmnXMMQhvjuAA0eeP7ddJklSZ10C1DXAIUmen+RJwInAJf1US7uqs87qf91KvdY6rXydtGvrNN1GklcDHwB2A86rqvfvbHsP8UmSxs13iK/TxWKr6lLg0i77kCRpLl5JQpI0kQxQkqSJtKxTvieZBW5dtgL7tR/w3ZWuxCpgOy2O7bQ4ttPirPZ2el5VrRtfuKwBajVLMjPXSTw9nu20OLbT4thOi7OrtpOH+CRJE8kAJUmaSAaoxTt3pSuwSthOi2M7LY7ttDi7ZDt5DkqSNJHsQUmSJpIBagFJfjXJN5I8mmRqbN1vtrMJ35jkF1eqjpPCGZbnluS8JHcluX5k2b5JLk/yrfZ+n5Ws4yRIcmCS6STfbP/n3t4ut61GJHlKkquTXNe203vb5c9PclX7//eJ9hqpq5oBamHXA78CXDm6MMmhNBfI/QfAscDvt7MMr0nte/8fwHHAocBJbRsJzqf5jIw6A7iiqg4Brmifr3UPA6dX1aHAkcBb28+QbfV4DwI/X1UvAQ4Djk1yJPA7wO9W1cHA3wFvXsE69sIAtYCquqGqbpxj1fHAx6vqwar6G+AmmlmG1ypnWJ5HVV0J3DO2+HjggvbxBcBrl7VSE6iqtlbVV9rH9wM30EyCaluNqMYD7dM92lsBPw98sl2+S7STAWrpFjWj8Bpiezwx+1fV1vbxncD+K1mZSZNkPfBS4Cpsqx0k2S3JtcBdwOXAt4F7q+rhdpNd4v+v09XMdxVJ/gx41hyr3lVVn13u+mhtqapKYjptK8lewKeAd1TVfdk+9zu21XZV9QhwWJK9gT8GfmqFqzQIAxRQVa9cwsucUfjxbI8nZluSA6pqa5IDaH4Jr3lJ9qAJTh+rqk+3i22reVTVvUmmgZ8F9k6ye9uL2iX+/zzEt3SXACcmeXKS5wOHAFevcJ1WkjMsPzGXAKe0j08B1nxPPU1X6SPADVV1zsgq22pEknVtz4kkTwVeRXO+bhp4XbvZLtFODtRdQJJfBv4bsA64F7i2qn6xXfcu4E002UfvqKrPrVhFJ8ATnWF5rUhyMXA0zRWntwFnAZ8BNgMH0Vzh/4SqGk+kWFOSHAV8Efg68Gi7+J0056Fsq1aSF9MkQexG08nYXFVnJ3kBTXLSvsBXgZOr6sGVq2l3BihJ0kTyEJ8kaSIZoCRJE8kAJUmaSAYoSdJEMkBJkiaSAUqSNJEMUJKkiWSAkiRNpP8PyKv2iEIfV3UAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 504x216 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfAAAACsCAYAAACEl/7eAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO2deZRtVX3nP7871K3hTTxAZBSjmAh2y1IcGxVnINqYTmLExIAxQbNi2644xrVsiTGGpFXsNEaDSh5qREknGJwnQtOmHYBunHAMPgR8gEwKyPCqavcf5xTcd+t8d9W5detWnVffz1q16t69zz77DL+zf+ec+/39dqSUMMYYY0yzaK31BhhjjDGmPnbgxhhjTAOxAzfGGGMaiB24McYY00DswI0xxpgGYgdujDHGNBA78AEi4riIuLbv+7cj4rjlLGvMahIRKSIeOsb+Di/77JTfPx0RpyxnWWOGISJOj4gPlZ8Pi4g7IqK91LIblUY78IjYGRF3lSd54e+sUfaRUjoqpXTxKNe53oiIUyPiS2u9HeMgIo6NiP8TET+LiFsi4l8j4jFrvV0rJSIujoi7B66Fj4+yj5TSCSmlc0e5zvVGU2/KR23XpT39fqZ+4YbtjoG/3xq2z0FSSj9OKW1KKc2Nap3rkYjYERFvGabt3nC3/NyU0hfWeiPM+icitgCfAP4QOB+YAJ4E3LMG29JehYHp5Sml9414nWadM0q7jogAokaTbSml2br9mNHQ6CfwHIOvVypeB26PiL+LiJ9ExK0R8TGxnp0R8Yzy81R5t3RrRFwJPGZg2YMi4h8j4qcR8aOIeEVf3WMj4ssRcVtE7IqIsyJioq8+RcTLIuIH5TLvKi+mqm1qR8QbIuLfIuL2iLg8Ig4t654YEZeWd+KXRsQT+9qdGhFXlW1+FBG/HREPB94DPKG8g75tiMPdFB4GkFI6L6U0l1K6K6X0uZTSN+C+4/q2iLipPE5/NGAz99lC+X3Qxv4hIq4vj/0lEXFUX92OiHh3RHwqIu4EnhoRvbK/H0fEDRHxnoiY6mvzmtJWfhIRvzfsTle9YYm+1/GlXb89Iq4ut/1L/dvR1+a+p7LBYwX86sCyWyPi/eX2XxcRb4nyVWhEPCQiLoqIm8v2fx8R2/ra7oyIV0fEN8rt+WhETGb27w8i4julXV8ZEY8qyx9ebvNtUfwU9h/72pxYLnt7uX2vjogZ4NPAQXH/E+VBQxzycbNSu744Iv48Iv4V+AXwQYobgLNiyLeaMfAEP2iDEXFURHw+ircFN0TEGyrWMThmPzgi/ld5zj4P7Dew/OOjeAtxW0R8Pfp++oyIF/fZyFUR8dK+uuMi4tqIeFVE3Fja7Isz+yZ9R2mLPyz368IF+4mCM8v1/zwivhkRj4iI04DfBl4bw7w1Syk19g/YCTxD1J0OfKjv++FAAjrl908CHwX2AbrAU8ry44Brq/oAzgD+N7AdOBT41sKyFDdDlwP/leIO+JeAq4Bnl/WPBh5P8dbjcOA7wCv7+kkUd9HbgMOAnwLHi317DfBN4Jcp7pYfCexbbtetwIvKfk4uv+8LzAA/B365XMeBwFHl51OBL631+RyDvWwBbgbOBU4A9hmofxnw3fLcbgf+ZcBm9rC3Chv7PWAz0APeCVzRV7cD+BnwH0pbmQTOBC4s+9oMfBz4i3L544EbgEeU5+7D5bY8VOzbxcDvi7pF57d/XcC7yvYHA23gieU+HD6w//f1sYxjdQHwt+W2PwD4GvDSsu6hwDPLPvYHLgHeOXDNfQ04qFz3d4CXiX37TeA6ipvpKNf9IIpr+ofAGyiux6cBt/fZ/y7gSeXnfYBHVV3/Tfhj5XZ9MfBj4CiKcaObs6eyzR62sZQ99tsgha3vAl5FcR1sBh43eE1V2N+XgXeUdvPk8nwuLHtweQxOpLi+nll+37+s/1XgIaWNPIXiRqX/nM8Cby73/cSyfh+xb8p3PA24CXhUuY3/A7ikrHs2hX/YVm7Dw4ED+8aGtwx17tfa+FZouDuBO4Db+v7+YNAQBo2BwnnNV50g8g78KvqcKnAa9zvwxwE/HljXnwB/J7b9lcAFfd8TcGzf9/OB14u23wNOqih/EfC1gbIvU1w8M+Xx+XVgamCZU9kADrzc14eXF8y15UV7IXBAWXcRfY4CeBY1HPhAP9vKtlvL7zuAD/TVB3An8JC+sicAPyo/nwOc0Vf3MJZ24L8YuBb+TJ3fhXVRDHZ3AY+sWOd910xfHwsOXB4r4ACK17dTffUnA/8itv15wP/r+74T+J2+738FvEe0/SzwXyrKnwRcD7T6ys4DTi8//xh4KbBloN1xNMyBj8CuLwbeXGFPy3Hgtw38PbyqPXs68JP7z/fAek+nwoFTPNTMAjN9y364b9nXAR+ssI1TRD8fW7Cb8pzfRd/NCHAj8PiKdjnf8X7gr/q+bwJ2l/vxNOD7FA9xrYF2OxjSge8Nr9Cfl1La1vf33mW0ORS4JaV0a82+DgKu6ft+dd/nB1G8ertt4Y/i7v8AgIh4WER8IopXrD8H3srAKyCKAWeBX1AYgNr+fxPbd/VA2dXAwSmlO4Hforgb3xURn4yIX1E7ureSUvpOSunUlNIhFE+3B1E8LUP+/GYpX1OeEcXPGj+ncEKw5znuX/f+wDRweZ+9fKYsH3ZbXjFwLbxxGW32o3gKqrKnHEtdC10KO1vYt7+leBInIg6IiI+Ur65/DnyI1bkWrkkpzQ9s48Hl51+neNK6unwt+wSx/kYwAru+pqJsOew3YHPfWUYbdc5yHATcWo5jCwza3G8OjL/HUjhcIuKEiPhK+Wr7Nopz329zN6c9f8tXNpfzHXuMvymlOyjeAhycUroIOIvibdeNEXF2FNqFFbE3OHDFnRQD5AIP7Pt8DbA9+n53Wya7KE7gAocNrPNHA8a8OaV0Yln/borXWEeklLZQOPc6YpF+rqF4HTTITygMuZ/DKF4xklL6bErpmRRG/V1g4WYnDbkdjSal9F2Ku99HlEW58wt5m3ohcBLwDGArxV037HmO+4/zTRR3/Uf12cvWlNLCoLHUttRhj+2OiP7tvgm4m2p7yrHUtXAPew7uW1JKC5qAt1Ici39XXgu/w+pcC4dGRP8Y138tXJpSOonipuJjFG+8YC+4Foawa1i83ys9DkuNv79Uc327gH1KncICgzb3wYHxdyaldEZE9IB/BN5G8VZiG/AphrO5nO/YY/wtt3Vf7re5v04pPRo4kuKN2mvKRYc+1nuzA78CeHIUsYRbKV5nA5BS2kUhVvmbiNgnIroR8eRlrPN84E/KNocA/7mv7mvA7RHxuihEQe1SpLAgdNtM8Rv0HeWT7x+uYN/eB/xZRBxRiiP+fUTsS2GUD4uIF0ZEJ4qQjiOBT5RPPSeVRnUPxU8PC08nNwCHRJ+obm8kIn6lFKocUn4/lOJ13lfKRc4HXhERh0TEPsDrB1ZxBfCC0l6OAX6jr24zxXG9mWLgemtuW8onw/cCZ0bEwpPpwRHx7L5tOTUijoyIaeBNw+01AF8HjoqIo6MQg50+sB3nAO+IQoTZjognlINeDnmsyuvrc8DbI2JLRLSiEK49pVxkM4X9/SwiDub+gWwY3ge8OiIeXV4LD42IBwFfpXiKem15vo4Dngt8JCImohBwbk0p7aa4LvuvhX3LMaMRjMCuq7iB+k62nyuA/xQR01GIJV/SV/cJ4MCIeGUUQs7NEfG43MpSSlcDlwF/Wp6/YynO5wIfAp4bEc8ubXgyCnHaIRQaiB6Frmg2Ik6g+BmhNkv4jvOAF5fXWY9iDPhqSmlnRDwmIh4XEV2Km5u72dPmhjrWe4MD/3jsGYd4AUBK6fMUQoNvUIgHPjHQ7kUUv098l+L3jlcuo68/pXhF8iOKAeqDCxWpCAl6DnB0WX8TxeCyMBC8muIp7XaKgfujtff0ft5BcVF+jmLweT/F7403l9vwKgpH8lrgOSmlmyjO9R9T3CXeQiHkWLiJuAj4NnB9RNy0gu1a79xOoVX4ahRK8K9QCBFfVda/l+J3s68D/xf4p4H2b6R42ruVwhY+3Ff3AQrbuA64kvsHzxyvoxBafaV8lfwFCmEiKaVPU7wCvahc5qJlrO+sgWvh8nJd36cQ6HwB+AEwGPP/agpR5KUUtvGXLD02LHWsfpdi4LyS4nj9T8rXmRTH7lEUor5PVrRdNimlfwD+nOJc3E7xNL09pXQvxQB/AsW1+DfA75ZPp1Bc/zvL4/4yCiXwwtPrecBV5avYJqjQV2rXVfx34DeiUFn/dWa52wZs7o/L8jOBeymc07nA3y80SCndTiEyey7FTyU/AJ66jG16Ybmft1Dc0H6gb53XULwBewOFo76G4sawVfb3Coox89ZyPRcuoz9Fpe9IRTjzGyme9ndRjBUvKNtsoTgPt1KMEzcD/62sez9wZGlvldFQiih/RDfGDBARh1PcjHWTY13NXoLteu9hb3gCN8YYYzYcduDGGGNMA/ErdGOMMaaB+AncGGOMaSArcuARcXxEfC+K3K/LCU0wZqzYRs16xzZqhmXoV+hRTEzwfYpwgGspQlBOTildqdpMRC9NMrOoPFriPkKVA7TrtUmt6ph9VQ6QKmehza1LrCizG6pNUpulls/dirXEORZtIqqXb7XmK8sB2qKPO75/w00ppf0rK1eZUdrokBugtqvW8rI8V6fsWiyfuw7UupJal7hu5DaRuQ5GVK7qZm+9hbk77hw2icyKqWujcgztioklO3rCyfmOGCs79cY3eb5zbUY07hV1anyrLpeXTGZ8a4kxUY177ahelyoH6Ii666+8TY6hK5lO9LHAD1NKVwFExEco4vDk4DjJDI+Lpy8qb01NVywNsVllT4SYqW4zv3nRJEoAzM1U5yiZndaHYHZTtWXunqq2pt3T1ZYxN6XHh9nq3WBWzL00N11tMHOT+kZsfqp61soQ5b2p3ZXl05N6dsJtU3dXll/89HcsOx3pKjAyG5W09OgV7eo6NdjGRLd6RV2dXyd6ok6sK01WLz8/JfoG5ier62anqvdvdkaUi+sG9LWze6a6fFbcY+2e0dfBbEXdT97+zoolx0otG1X22dn/gRVLw/wB22XH9+5bPVbes63aPu/dLMa9TXp82y3O01x118xOifFtWju+NFld15qqjpCb6FWXT/XulX1sEnVbetXj3pZudfm2ibtkH9s6v6gs/8uj/0mOoSt5hX4we+bPvZb78wwbsx6wjZr1jm3UDM1KnsCXRRTznZ4GMIl43DRmDbGNmvWM7dMoVvIEfh17Jsg/pCzbg5TS2SmlY1JKx3RZKr2yMSPFNmrWO0vaqO3TKFbiwC8FjoiIB5eTYLyAleWXNWbU2EbNesc2aoZm6FfoKaXZiHg5RZL8NnBOSunbuTbRalUK1mKqWrEVXS2sSUIIlOqq0Ns5ZWxNtblYVVYZq7pX66q7nsy6lNpcqtDXTKs7HMPYqEQpyodQb2u1uYqsqK9Cl6ryIZTuWhms1lVzPbm6uuuqex2ssU3XtdHodioFa2nb5srl5zLixLmJanubFyr0eRmVI7vIqNCrx5i6Y2tRV09tXnfcA61Cr7t8KzNzaN0+YIW/gaeUPkUxhaUx6xLbqFnv2EbNsDgTmzHGGNNA7MCNMcaYBmIHbowxxjQQO3BjjDGmgax6Ipc9aLUq06MqtblK+wiQekKF3quWSqZuvby/kFNQ1swVPJRKUzSQiuCMgrGmYncYlabK49tIKtKjSrW5SJcKmVSqqo3oQ+ZOz7SRcwKIqItsLvSabWpHb4DO8V9XbZ5Tuu8NjyudTmV6VKU2n53RQ/zcpFChC+F6Equaz3iR2vnTpTo9M77JsXJ0UTa6jSoX6V0z42QuT7pibzBpY4wxZsNhB26MMcY0EDtwY4wxpoHYgRtjjDENxA7cGGOMaSB24MYYY0wDGW8YWbtFzCyezEROTCJCxYq66liH+W51fIJM3N/VsQMyqb8Mp1ATpsguMmEWIpxClGdvxUQIRojylihvt3SYQy7ErFFEVId5qRCvTBiZDDFT5R1hWKo824eaGEWFl2kDGlXYZC5UrW64Ue3wJMiHWjaE+U6Le/edWlSuxjcVKgYwO1V9PmYnq8vnJsR4qOdLkaFn8rzK8DLdhxzHVIjXEONbR9R1Yq56XTXDziA/0YluY4wxxpjGYQdujDHGNBA7cGOMMaaB2IEbY4wxDcQO3BhjjGkgY5/MZH7zYgVlUspYMTEJZNTmPaU2FyrNjApd1UllZU3FZa6ufrlWMIaoawllpVJKKmUlaJVm04gIoioqQkwoMpQKXajKQ1wHdHQfSj2eRN+pU2+SE4B5NZmJKNfLyy5qt6mrTi/aVNhvw6InUie4Z9ti+5ERMxmFuFSbq/Je9Xpyfcx3q4/vfEdE2YhysuObGMdEuRqrcip0NdGIVqdXl3eFaj3XRw4/gRtjjDENxA7cGGOMaSB24MYYY0wDsQM3xhhjGogduDHGGNNAVqRCj4idwO3AHDCbUjomt3xqBXMzE5XllcsL5TjkcpurnMD18vsW6xJ9C3W6zJGeU6ErNWbNHOm5PM9KpdkW6+q0q5WSXVEOMNGalXVrSV0bJYKYqDjxoXKLa/uprTav6hdIGRU6Yh4BlNq8Zjnoa0rNI6AU4ur6KPoX5TWjLnLRGJVK5szpGxd1bDS14N7Ni8+HGmPUcQU99km1+eKhuyzXx1yNoapcqtA7WqHdEue8rVToorw7VC50ta7qsTKXCz2nUJfbVbvFYp6aUrppBOsxZrWwjZr1jm3U1Mav0I0xxpgGslIHnoDPRcTlEXFa1QIRcVpEXBYRl+2evXOF3RlTm1o2em+6a8ybZ0zeRvvtc/Zuj6Hmflb6Cv3YlNJ1EfEA4PMR8d2U0iX9C6SUzgbOBtiy6eBmpT0yewO1bHRrZ3/bqBk3WRvtt8+Z/Q61fZr7WNETeErpuvL/jcAFwGNHsVHGjArbqFnv2EbNsAz9BB4RM0ArpXR7+flZwJtzbVIrmJ1e3KXKqZxEfl/QClilEJeKy4wKfU6qLtU21SuHjHJd5BDWKs36udCVSrMryieEshJgIqNQXyuGsVEioFtxgoXaPESOdECq0FVuc6k2V0pzIIk5AeYnRLlaPnOtqetQ5t+W5bKL2hEccvnciJa5RtaKujaa2rB70+Ljq6JWsse85nil1OZqnMy1SV2h+BZq85ZaHmh3qscepTYfZZRNr11d3hGK8snWbtnHuFXoBwAXlANYB/hwSukzK1ifMaPGNmrWO7ZRMzRDO/CU0lXAI0e4LcaMFNuoWe/YRs1KcBiZMcYY00DswI0xxpgGYgdujDHGNBA7cGOMMaaBjCIX+rJJbZjdtDguRE5mkrm9UKEqcqIRNTFJJgRirifWVTeMLJPsX07iUDNcLDLJ/nW4WHXYggoJa+JkJrWJIHoVJ1iFi+UmM2mLcLG2MGwRLqZCxQDSRL02cmISMTkQZEIw5WQ/1eUq9BMy14Gc9EItr6+DqKpr2CNMasHumepytbxcV80Q1qHGN1Un+ogJNfFSZqIRMfZNiPCynhrfhgiTVZOZ9MR4mAsVGyaMrGHma4wxxhiwAzfGGGMaiR24McYY00DswI0xxpgGYgdujDHGNJDxqtBbwe6pxfcMWkGZmcxEqrery6U6PZuIX61LLS/Um7kJBVSyf6neFCrNjAq92xWqS6XS7FQrKCfbOhH/VKauUbQCJhafYGmLrcw9sFCbJ6FOp6MU4hkVulCbz/VUeXUfqhy0elyW14zSgIyiXV1TNdXNoK6R9TfBSY7UgrmpqnIxjmjTkXXzItJFnb+c8l+rzavHHjWOdSd0lMuEGK9UuZqAJDe+qSibqfa91X2oyU8yk5nk6hR+AjfGGGMaiB24McYY00DswI0xxpgGYgdujDHGNBA7cGOMMaaBjFmFDrunK5SrQuA7TC50md9X5v3N9FFTba7yqs/3tNJ1XuT+RagxW0Lx2RFKc9BqzMmaKs3pTrXiEvYiFXoEaXL5udBTO5MLXSjUk1Cbq/L5TC50ldtcqcrnRX5/le88V6fsXS0/TMSHvgarrwOlboZqJXMI9fZ6JQXMTi3eZjlWZvZPqdDVPAx6foZM/nmV27ym2lzlNQfodWuObzLKRivd1fim1OaTQlE+GbqPydDjq8JP4MYYY0wDsQM3xhhjGogduDHGGNNA7MCNMcaYBmIHbowxxjSQJVXoEXEO8BzgxpTSI8qy7cBHgcOBncDzU0q3LtlbC+amFitUdS50vSqpoFT5fVV5ToVeM7e5UpvLvOaQyRVcT23eFUpM0CrNqU61UlKpzXNK8+lWfQXlqBiljaZWMD9VYRRKhZ7J159kLnQVQSGU4yLiAmB+op4KXSvKMyr0nioXanO1fE6FLq4dmfNcqZszebknJxbbbyvGo0IfmY22YG66Yh9lJE9m/9T42q6nNleRMQDt9mjU5moMA602V+ObynmeG99U3bTIha7Gw1y+c6Vcz7GcJ/AdwPEDZa8HvphSOgL4YvndmLViB7ZRs77ZgW3UjJglHXhK6RLgloHik4Bzy8/nAs8b8XYZs2xso2a9Yxs1q8Gwv4EfkFLaVX6+HjhgRNtjzKiwjZr1jm3UrIgVi9hSSonMpLoRcVpEXBYRl83ededKuzOmNnVsdPfsL8a4ZcYU5Gy03z7n7rhjzFtm1jPDOvAbIuJAgPL/jWrBlNLZKaVjUkrHdKZmhuzOmNoMZaPdzvTYNtBseJZlo/322d60aawbaNY3wzrwC4FTys+nAP88ms0xZmTYRs16xzZqVsRywsjOA44D9ouIa4E3AWcA50fES4Crgecvp7PUgtmKB5w0xGQmMvSsZoJ+FRJWtKkul+FlamISESoGEL3qsImOCrMQ5ZOZMIvprgibEGEWMyKMbFP7HtlHrm61GaWN0grmJxeHkSkbJTOZiQoxmxdt1MQkKRNGJicO6dYLF1OhYkVdvTbDTOozJ66p1BOhS+K66U1mQh3XNoxsNDbaSqTJqjAysR+ZMVRN5BIi9KslwsvamYlGOiL0TE2wVHdiEqgfDqvGt5mOHsNUWFjdcLGZlu5jMuqHkS3pwFNKJ4uqp9fuzZhVwDZq1ju2UbMaOBObMcYY00DswI0xxpgGYgdujDHGNBA7cGOMMaaBLCliGyUpYHayomIYFbpQREoVuljXfEYhrpTrcnISobhUE5NATm1erexUanOlNAet0tzUrVZEbu7cXb18RqW5uV3dpmmkCGanKoxIKMqzk5nISIl6yvHcZCZzYjIerUIX6xFK86KuZvlk9fUxl1Ghz1cpq4GYFJP39MR10NOT6mzpLbbR9phU6COjlWhNLd53MdeOVqejFfgtoUJXE5N0RDnoyUmUCr3uxCQwhNpcTECSnayp5qQlSm2eU5rnFOoKP4EbY4wxDcQO3BhjjGkgduDGGGNMA7EDN8YYYxqIHbgxxhjTQMaqQqcFc9OLlY96nke9qtoqdLW8UJoDoOq6IlewKO90da7gurnNZyZEnvIJrWDc0q1WiKtyldd8a/su2cfW9t4xVWxqw+xMhRHJSIn6KnSVC13Z7lAqdKU2V7nTc7nQ1bqU2lyUz09lrjWhQu8otflk9XWwdVJHQ2zvLZ4qttPSCur1SARMVByTEIpyVQ7QErnQ1TFRavNOW49vPVHXa4tc6EqdnlGI11Wbq5znuUiaTaJuWijHVfmoc6H7CdwYY4xpIHbgxhhjTAOxAzfGGGMaiB24McYY00DswI0xxpgGMt5c6K1qhWpSIluhkizqRB9CbS7XlVGhh8ht3hblSm3eFYpyqJ/bXKnNVV5zgM1Cba5ynm/tVKvNt7YXq3gX2LKX5EKnFcxOLTYuZaPZfP0qf7pUm9crh0z+9Nrq9EwfIoe5ym2u1OZJ5DUH6EyKHNhT1UriLZPV9r5PhdJ8gf16dyzuN/Q2rUdarXmmKvK9K7V5JkiCtlCbq/KuKs+o0LstMaeDUJVPCnV6Lk+5UpWrNkptrvKa5+qUqlyVK3U6wIxV6MYYY8zGwA7cGGOMaSB24MYYY0wDsQM3xhhjGogduDHGGNNAlnTgEXFORNwYEd/qKzs9Iq6LiCvKvxNXdzON0dhGzXrHNmpWg+WEke0AzgI+MFB+ZkrpbbV6ayXmpyrCClSoQyYEQoWFhQgjC5GIXy0P0BZtuiJcbEIk4u9lwshUuNhUp7pcTUCiQsUAtnVVWFi9cLFtmTCybS1dNwZ2MCIbTS3YPb3Y8HSoY35dVcjJTIYII1Nt6k5yokLFAOYmRLiYmIBETkwiQsUApqbFBDpT1Xa9fbJ68pwHVISKLXDgxM8WlXXHF0a2gxHYaCsSmyrCyFpDTGaiJi1phzh/YvmJlj6vEyLETLVRoV+5MDIV4jUtJjPRE5PoMLItreqxUk5aEiLsTJQX69LHUbHkE3hK6RLgltprNmZM2EbNesc2alaDlfwG/vKI+Eb5amifkW2RMaPDNmrWO7ZRMzTDOvB3Aw8BjgZ2AW9XC0bEaRFxWURcNnf73jFntGkEQ9no7F22UTM2lmWje9jnz6pf5ZqNyVAOPKV0Q0ppLqU0D7wXeGxm2bNTSseklI5pb54ZdjuNqcWwNtqZso2a8bBcG93DPrdOjXcjzbpmKAceEQf2ff014FtqWWPWAtuoWe/YRs1KWVKFHhHnAccB+0XEtcCbgOMi4mggATuBly6rtxZEhQpdKiUzKvQQKvSWUlYKtblSmoNO0t/tiAT9QoWuFOW5OjU5iVShi4lJQKvNt3eqFbuqfFtbv17eJlSa42CUNppasHumwvDUZCYZG1WTltQtnxflAEmozaVyXKrTtVo59UQEh5icpNOrNzEJ5NTm1dEND5y6vbq8t1hpvsAhE4s1ZGqyjVEzKhtttxJbesufOEip00FP5KLU5h2hTu+JCUhybaaEQrwnlNhKUQ4ZFbosrzcBSa7NlpZStIvJTCIzsVXmXCmWdOAppZMrit9fuydjVgnbqFnv2EbNauBMbMYYY0wDsQM3xhhjGogduDHGGNNA7MCNMcaYBrKcXOgjIyLRm1qsulYq9Fwe35ZSoYs2HaUoz6jQVR7fnsp5LtSY0x2toJwRdUpVvqktckZ3tApc5TZXavN92wrzAaEAAAW6SURBVEKFnlGab21ppX2TSC2YrQgFl2rznApd3B5rFbpQjg+RC12pyucnhL2rcqDVE9eOUJtPT1bb9JZJrfJVuc2V2vyg3m2V5YdN3CT7OLR786KyiYwqeD3SjvnKSBQ17rWECrxYlxgrRRul2FdqdtCqclU+KcaRXJ7ynmijVOXDqNBVDnO1rs2htkmfj+nITf5RjZ/AjTHGmAZiB26MMcY0EDtwY4wxpoHYgRtjjDENxA7cGGOMaSBjVaG3WvNMVyhRW0J8l1Oht1XOc9FG5TWfyORCVm0m20IpKRTlU2J50KryTR2hNhd5zZXSHGCbqFO5zZXafHtGCbq9nUnY3SCKXOgVNjSMCl2pzaU6fQgVeleoWrsismOi2qbbaj1Ab1LYe6/aHrZOVkdQ7NPTNvqAXnXkg8ptrtTmh3UX5ztf4PCKqIse48mFPiraMc+2icXXZ4v6KnSlXO8KVblaXinHc+tSyvFJERWglgetHp+USvB66nTIqNDF9iq1eU5pPh1iooIMfgI3xhhjGogduDHGGNNA7MCNMcaYBmIHbowxxjQQO3BjjDGmgdiBG2OMMQ1krGFk7VZi21R1iEkVKqk+6BCzjpDvT4jk+WrCklwbFRamynOJ+FUY2eZ29XHaKkK/tojlAba1RBiZCBdTE5PkQsW2tqZkXZMoJjNZfhiZCgkDQEy4o8LFUOUdHU4ZIvyr3RETUkyISSQmdJjOtKjb0qu2ue0iXGw/ESoGcOBEdbjYIRPVYWFVE5NAdajYAod1Ni0qm4hb5fLrkU7Ms62z+PiqEK92LoxMhJ6pNiokTJXn6nQYWfVYmQtVqxsuJpcX5QDTavIVcdxVuFguVGy6NSHrFH4CN8YYYxqIHbgxxhjTQOzAjTHGmAZiB26MMcY0EDtwY4wxpoFESlrhOvLOIn4KXF1+3Q+onpFg9dmofa9F/w9KKe0/xv5WhG10w/Vt+xyejWQna9m3tNGxOvA9Oo64LKV0jPveWP03iY1qJxu176ax1sdqo9rJWh/3fvwK3RhjjGkgduDGGGNMA1lLB362+96Q/TeJjWonG7XvprHWx2qj2slaH/f7WLPfwI0xxhgzPH6FbowxxjSQNXHgEXF8RHwvIn4YEa8fc987I+KbEXFFRFy2yn2dExE3RsS3+sq2R8TnI+IH5f99xtj36RFxXbnvV0TEiavRd9NZS/ss+7eN2kazbJQxtOzPNioYuwOPiDbwLuAE4Ejg5Ig4csyb8dSU0tFjCAXYARw/UPZ64IsppSOAL5bfx9U3wJnlvh+dUvrUKvXdWNaJfYJt1DYqWCc2Oi77BNuoZC2ewB8L/DCldFVK6V7gI8BJa7Adq05K6RJgcD7Ek4Bzy8/nAs8bY99maTaMfYJttKHYRm2jwNo48IOBa/q+X1uWjYsEfC4iLo+I08bY7wIHpJR2lZ+vBw4Yc/8vj4hvlK+GVuW1U8NZa/sE26htNM9a2+ha2yfYRoGNKWI7NqX0KIrXT38UEU9eqw1JRQjAOMMA3g08BDga2AW8fYx9m+VjG7WNrmfWjX3CxrbRtXDg1wGH9n0/pCwbCyml68r/NwIXULyOGic3RMSBAOX/G8fVcUrphpTSXEppHngv49/3JrCm9gm2Udvokmz0MRRso8DaOPBLgSMi4sERMQG8ALhwHB1HxExEbF74DDwL+Fa+1ci5EDil/HwK8M/j6njB4Et+jfHvexNYM/sE26htdFls9DEUbKMAdMbdYUppNiJeDnwWaAPnpJS+PabuDwAuiAgo9v3DKaXPrFZnEXEecBywX0RcC7wJOAM4PyJeQjGr0PPH2PdxEXE0xeumncBLV6PvJrPG9gm2UdvoEmykMRRso9ntcyY2Y4wxpnlsRBGbMcYY03jswI0xxpgGYgdujDHGNBA7cGOMMaaB2IEbY4wxDcQO3BhjjGkgduDGGGNMA7EDN8YYYxrI/wfsdIQCMASDMQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 504x216 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "n = 20  # nb samples\n",
    "xs = np.zeros((n, 2))\n",
    "xs[:, 0] = np.arange(n) + 1\n",
    "xs[:, 1] = (np.arange(n) + 1) * -0.001  # to make it strictly convex...\n",
    "\n",
    "xt = np.zeros((n, 2))\n",
    "xt[:, 1] = np.arange(n) + 1\n",
    "\n",
    "a, b = ot.unif(n), ot.unif(n)  # uniform distribution on samples\n",
    "\n",
    "# loss matrix\n",
    "M1 = ot.dist(xs, xt, metric='euclidean')\n",
    "M1 /= M1.max()\n",
    "\n",
    "# loss matrix\n",
    "M2 = ot.dist(xs, xt, metric='sqeuclidean')\n",
    "M2 /= M2.max()\n",
    "\n",
    "# loss matrix\n",
    "Mp = np.sqrt(ot.dist(xs, xt, metric='euclidean'))\n",
    "Mp /= Mp.max()\n",
    "\n",
    "# Data\n",
    "pl.figure(1, figsize=(7, 3))\n",
    "pl.clf()\n",
    "pl.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')\n",
    "pl.plot(xt[:, 0], xt[:, 1], 'xr', label='Target samples')\n",
    "pl.axis('equal')\n",
    "pl.title('Source and target distributions')\n",
    "\n",
    "\n",
    "# Cost matrices\n",
    "pl.figure(2, figsize=(7, 3))\n",
    "\n",
    "pl.subplot(1, 3, 1)\n",
    "pl.imshow(M1, interpolation='nearest')\n",
    "pl.title('Euclidean cost')\n",
    "\n",
    "pl.subplot(1, 3, 2)\n",
    "pl.imshow(M2, interpolation='nearest')\n",
    "pl.title('Squared Euclidean cost')\n",
    "\n",
    "pl.subplot(1, 3, 3)\n",
    "pl.imshow(Mp, interpolation='nearest')\n",
    "pl.title('Sqrt Euclidean cost')\n",
    "pl.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Dataset 1 : Plot OT Matrices\n",
    "----------------------------\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAe4AAADQCAYAAADBEII/AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOydd3gU19XG36te6Agwokl0JLApEh2bARyKA6Z5AVNEcGy8TrETQgiJRWzJDrFJIU78reWSIMAGRoBlbINN8WCqJCRM75gmECCJJprq/f44O+zsaiXtrna1u9J9n0ePpNnZ2Zm975zf3DtnzmWccwgJCQkJCQl5h3zcvQNCQkJCQkJCtkuAW0hISEhIyIskwC0kJCQkJORFEuAWEhISEhLyIglwCwkJCQkJeZEEuIWEhISEhLxIAtxuEGNsGWPsLePfQxhjJ21ZV0jIVWKMvcEYW1nDnynOAyGXizHGGWMdjX9/wBiLt2VdT1atAzdjbDZj7DBj7D5j7CpjzMAYa2R87QPG2F3jTxFjrFjz/yYr2xrKGCvTrKP+DHDW/nLOd3LOuzhre0L2y5meqY0S54H3yFO9bPRQdhXrLDPul9ZjB525H5zzlznnic7cpjtUq8DNGJsH4B0A8wE0BNAfQDsAWxhjAcZGq8c5rwfgLwDWqP9zzkdXsNkrmnXUn701ckBCLpeLPONWMcb8XLBZcR54uDzVy3b68V0Ljz3hqv3yZtUacDPGGgB4E8CvOOffcM6LOefnAegARACY4YLPPM8YG6H532y4kTE2mDG2hzF2izF2iTE228o2zK5EGWO9GGP7GWMFjLE1AIIs1v8pY+yAcZt7GGOPa177A2PsrPG9xxhjEzSvzWaM7WKM/Y0xdpMxdo4x5pHgqSm5wjOMsTDG2FfG9rnBGNvJGPMxvmbWtoyx1Zqh4tmMsV0W29IO8T3DGPuBMXbH6KU3NOtFGNd9gTF2EcB3xuX9Nf47yBgbqnlPJGPse+O+bAEQZu+xarYlzgM3y8O8PJQxls0YW8AYuwpgFYBNAMI1PelwO/elXI9d6zvGmC9j7I+ads9ijLWxsh2zWy6MsfmMsRzG2BXG2ByLdQONPrnIGLvGaMQi2PhaY+N3k2v00VeMsdaa925njCUyxnYb92czY8zhc8xStQbcAAaCTu712oWc87sANgJ4uiZ3hjHWDmTWfwNoBqAngANVvCcAQCqAFQCaAEgBMEnzei8A/wUwF0BTAEkANjDGAo2rnAUwBHS1/SaAlYyxlpqP6AfgJChIvwvgE8YYq9aBerdc4Zl5ALJBbd4CwB8B8Kra1gbdAzALQCMAzwDQM8bGW6zzFIBuAEYyxloB+BrAW8bP+x2AdYyxZsZ1PwOQBfJCIoA4+w7TNonzoMbkaV5+zPhaO5BvR8N81OaKA/tTmX4LYBqAMQAaAJgD4H5lb2CMjQKdF08D6ARghMUqfwXQGeTZjgBaAVhkfM0HwP9Ax9cWwAMA/7F4//MAfgagOYAA42c5RbUJ3GEA8jjnJVZey4HjPYpw4xWn9ifUhvc9D2Ar53yV8eo3n3NeacACDW35A1hqfM9aAPs0r78EIIlzns45L+WcJwMoNL4PnPMUzvkVznkZ53wNgNMA+mref4Fz/hHnvBRAMoCWoBOyrsoVnikGfa/tjG24k9OEAFW1baXinG/nnB82tu0hUC/mKYvV3uCc3+OcPwD1sDZyzjca37MFQCaAMYyxtgBiAcRzzgs55zsAfFnFLojzwLPlaV4uA/Bno78e2PGZv7PwWLKN7/s5gNc55yc56SDnPL+K9+gA/I9zfoRzfg/AG+oLxgu5lwD8hnN+g3NeALq9MBUAjD5exzm/b3ztbZQ/H//HOT9lPH4ZdAHgFNUmcOcBCGPW76e0NL7uiK5wzhtZ/Nyz4X1tQFf+9igcwGXjyaHqgubvdgDmaY1t/JxwAGCMzdIMH94C0B3mJ+xV9Q/OuXo1Ws/OfaxNcoVnlgA4A2AzY+xHxtgfjMurattKxRjrxxhTjENztwG8jPLB+JLm73YAnrPwymDQcYUDuGnh46r2RZwHni1P83Iu5/yhA5/5NwuP2ToS5KjPtOeM9hiaAQgBkKXx0TfG5WCMhTDGkhhjFxhjdwDsANCIMear2cZVzd/34USP1SZw7wVddU/ULmSM1QMN02xzwWfeAzWuqsc0f18C0MHO7eUAaGUxbNfWYptvWxg7hHO+yjgk+RGAXwJoyjlvBOAIAG8fAnSlnO4ZznkB53we57w9gHEAfssYG46q29bMS4wxrZcAGtreAKAN57whgA9Qvm21gfQSgBUWXgnlnP/VuC+NLXrMbeG4xHngfnmSlwFzL1r7315Znh++MELUKEd9pr0Prj2GPNDwd7TGYw05JfYBdBuhC4B+nPMGAJ5Ud83OfXBItQbcnPPboPtZ/2aMjWKM+TPGIkBDFNmg+zHO1gEAU42fFQNgsua1TwGMYIzpGGN+jLGmjLGqhkr2AigB8GvjNifCfIjvIwAvG3tfjDEWyihpqT6AUNDJkQsAjLGfgXoaQhXIFZ5hlDTV0RjUbgMoBQ0bVtW2BwFEM8Z6MsaCoBm2M6o+gBuc84eMsb6gIejKtBLAWMbYSEaJO0GMEnxac84vgIbN32SMBTDGBgMYa++xaiTOAzfLw7xsTdcANGWMNbR3P4w6BSDI2M7+AF4HEKh5/WMAiYyxTkZPPM4Ya1rFNmUAsxljUYyxEAB/Vl/gnJeBfPZPxlhzAGCMtWKMjTSuUh8E9luMsSba99aEag24AYBz/i4ogeJvAO4ASAddiQ3nnBc6uFltJqT6oyZixIOu8m6CTprPNPtyEZQoMQ/ADVBwq/TRBs55EeiKebbxPVOgSTbhnGcCeBGUBHETNIw12/jaMQB/B51U1wD0ALDbwWOuM3KBZzoB2ArgLqgt/o9zrtjQtqcAJBjfexrALvPN4hUACYyxAlCCjFzFcV0C8Kzx2HKNxzQfpnP+eVCS1g1Q0FlexXGJ88DD5SlermDfToDyMn40Dj1XlFX+ewuP5Rnffxt0DnwM4DKoB67NMv8H6JzYDDr2TwAEV7FPmwAsBT2Fccb4W6sFxuVpxuHwraBeNozvCwb1zNNAw+g1JmZ+m0JISKimxBhbBiCbc/66u/dFSKg6El6uWdWqHreQkJCQkFBtlwC3kJCQkJCQF0kMlQsJCQkJCXmRRI9bSEhISEjIi+SKyQgqVFhYGI+IiKjJjxTyAGVlZeVxzptVvab7JLxZNyW8KeSpqsybNQruiIgIZGZm1uRHCnmAGGM2Vwhzl4Q366aEN4U8VZV5UwyVCwkJCQkJeZEEuIWEhISEhLxIAtxCQkJCQkJeJAFuISEhISEhL5IAt5CQkJCQkBdJgFtISEhISMiLJMAtJCQkJCTkRRLgFhISEhIS8iIJcAsJCQkJCXmRBLiFhISEhIS8SALcQkJCQkJCXiQBbiEhISEhIS+SALeQkJCQkJAXSYBbSEhISEjIiyTALSQkJCQk5EUS4BYSEhISEvIiCXALCQkJCQl5kQS4hYSEhISEvEgC3EJCQkJCQl6kKsHNGGvDGFMYY8cYY0cZY68alzdhjG1hjJ02/m5s1ye/+y6gKObLFIWWCwnZIOFNIU+VS7wpfClklC097hIA8zjnUQD6A/gFYywKwB8AbOOcdwKwzfi/7YqNBXQ6lG41GlFRAJ2OlgsJ2SaXevNRkBTeFLJfzvemxpclJRC+rMvinNv1A+ALAE8DOAmgpXFZSwAnq3pvnz59uFYFG77j90LD+NW58ZyHhXH+3XdcqPYJQCa302eO/DjTmzfXf8fvh4bx/F8Ib9ZmeZ03v/uOlzQJ43ufjuclTYQva7Mq86Zd97gZYxEAegFIB9CCc55jfOkqgBYVvOclxlgmYywzNzfX7LWAkRJOD9ejRVIirk3SA5Jkz+4ICT2Ss71Z+qSEjD56NHk/EdALbwo5Lqd6U5JQ8nM9+m9JRFpPPbI7CV/WRdkMbsZYPQDrALzGOb+jfc14dcCtvY9z/iHnPIZzHtOsWTOz1wJ2K3h8jwEHx8Wj3goDfvxEsbYJIaFK5QpvBqcpiMk04PIL8YDBUP7eopCQDXK6NxUFgf81oPD38eidbsCONxVkZ7vwAIQ8UjaBmzHmDzLfp5zz9cbF1xhjLY2vtwRw3a5PNt6fYbKMbikJ2PUrGS1+rRPwFrJLrvJm8Gwd1j4n4/zPEgBZNr/nLSRkg5zuTfWetiwj8J0ElK2WMWGVTsC7DsqWrHIG4BMAxznn/9C8tAFAnPHvONA9HNu1ZAmwcCEgSQgIAKQECcfGLUTpO0tw5IhdWxKqo3KpN/+wEOcjJRQVgYbJFy6k5UJCNsgl3tTETAAI/akE/HEhBuxegpUrIeBdh2RLj3sQgJkAhjHGDhh/xgD4K4CnGWOnAYww/m+75s8HFi9GXooCzmnYPGbrYpydMB/r10PAW8gWuc6bf12MjpcU3LsH6uksXkzLhYRsk/O9aYyZUBRcvgwaGfrnYjR7Zz5CQiDgXYfkV9UKnPNdAFgFLw93+JMlCbc+lBE8Q4dTa/Xo/J0BTJYxbJCEq58B640DS927O/wJQrVcrvTmDYOM8XE6nMnTA7sNNFwuEtSEbJRLvClJgCyjZJIOZ3ro0fyQAf7rZdSTJMTdBpKTCd4zZgCtWzu860JeILdWTms4XkLuJD26yIk4NUwPPpSGzZ9/HmjbFqLnLeQ2BY6SkBmjxxNfiKxyIQ+SJMHnFT2e2pGIPY/rsTeIfNmwIRAXB9HzriNyK7jZdgXtNhlwfmY8Wn9lwJ63jcPmAt5CblZoBmWVp48UWeVCHiRFgU+SAWV/ikf/Hww4laRg7156yRLely+7d1eFXCf3gXvuXGDCBDBZRrvkBJxKlNHnrQk4P2qugLeQezV3LtjECVg/RcaO4cas8gkTyLNCQu6SMWZCluHzVgL81suYljIB/r+caxXeK1YIeNdWuXeSEaoeBMaAnj0BXz+O/Hzgyy8h4C3kXnEOPz9QaUnj/0JCbpfGh76+gL8fR+PGwObNEPCuQ3IfuJOSwD9PRfEEHYr+sAhsig5+G1JRsCQJP/wg4C3kRiUl4cGqVDz7mQ6Dvl1Ez86mpgJJSe7eM6G6rKQkIDUVXKfDhbhF4DodWGoqIr5NQrduAt51SW7tcec/LiGtlx4B7ySi6AU92DAJQ4cCTz4JAW8ht+purIR9ffR48nuRnCbkQZIk5Ov0aLc8EUeG6FEyRIKvLzBpEgS865DcCu6wwwoGHjJg59B4lL5vwIONChiDgLeQ29XkICWnff9kPLhIThPyFCkKwmQDrrwQj/bfGrD9zzRTmIB33ZL7wG0s3+e7VkarTxKwTicDU3Rm8B4yRMBbyA1SFPhO0yF1moztwxJw8wNR8lTIA6QpeRr+cQKy/y5jwFJdhfBOS6O3CXjXPrkP3Jryfe3bAwP/JGHXkIW4/vsluHePEtYkqWp4Hz3qtiMQqq0yejP/cRoevx4tSp4KeYAsSp52eVlC/ksL0VZeAllGOXh/+62Ad22V+8BtLN/Hv6NeTPsLCoalL8bOfvOxfDlsgnebNsC6dQLeQk6W0ZtdrpA3izeLkqdCHiBNyVPOASgK2q5cDD5vPk6fBlJSbIN3cLCAt7fLfeCWJBSukFH4rA63Xl30aNh84J8k3LgBm+A9fbqAt5ALJEm4/ZGMIf/RYeh3i9A5XidKngq5X8aSp2XP6XBg3CKUPUe+7PKyhDFjgFOnbIP37NkC3t4utyanFQ+WcHiwHo3eS8St5/WPhs2nTUOl8P7qKwFvIdfKdwRllT+1IxFHh4isciEPkSTh/iw9en2ViIzeetzsSb6MjYWAdx2SW8Fdb5+CPvsM2DcqHgEfG3B1lXHY3ALe9++bw3v/fgFvIdeq3j4FsVmUVd51u8gqF/IQKQrqrTCg4LV4PL7HgK1/UnDzJr1UGby7dhXwrk1ye1a5T4qMrnICNv9cRoOf66zCOzlZwFuoBmX05t7XKKt8/RSRVS7kAdJkldf/ZwIeLJPxzHKdTfCePLlqeF+54rYjE7JTHpFVXr8+MPwtCVkjF+LuG0tw8SKtIuAt5BYZvRkwkoYhz7YVWeVCHiCLrPKmkyUU/24hen+3BMnJcBjeasLa8uUC3t4it2eV3/+aejH1MxUM2rEYh0fOx8qVsAnegwcLeAu5QEZvtr9A3ow4p4CLrHIhd0uTVV5QAEBR0PD9xWj01nwUFsJueKen07qNGgl4e5vcmlV++yMqunJx9qJHw+Yj3pbQoAEqhLf2nvewYQLeQi6QJCH3fRmN5lJW+eQUHS7/Q2SVC7lZxqzy0kmUVV4yiYbNm06WMGsWqoR3aak5vL/5xjq8xbC558utyWn1x0m4MFqPtsmJuPiM/tGweVwcKoR3fr5t8NY+533smDuPUsgb1XiihB/6UVZ5ZoweZ9oIaAt5gCQJ/GU9hmxPxJ7H9Tj+GPmyZUtUCW9Ztg3eQUEC3p4ut4Lb53sFXbcbcHxyPJqmGHDoX8Zhcw28P/3UHN5Tp9oG78BAE7zXrhXwFrJPfjsVxBprlcdkGqgIi5CQu6Uo8PvIgJKF8YjNMmDfuwqOH6eXKoP36NEC3rVJVYKbMfZfxth1xtgRzbI3GGOXGWMHjD9j7P5k46TwTJbRZU0Csn4vo/OCCcgZNxeACd7165vDu0OHquH99dcC3nVBrvbm5X9SVvna52Q8+c8JtFxIyAa5xJtGX0KW4feXBPiulTFlzQTcnzm3Snj37SvgXZtkS497GYBRVpb/k3Pe0/iz0aFPN04K7+ND0PXx5bh8Bdi1i152FN5ZWeXh3bq1gHct1DK40JutW5v+LTN6VUjIRi2DK7yp8WFAABDgzxEaSrHtxAlabgu8bb3nLeDtmaoS3JzzHQBuOP2TjZPCl0zSofRPi+AzVQe/Dam4+MckbNvmfHhPny7gXdvkSm+WrkuFzzQdhn9PyWny1FQU/yfJ6R8lVDvlEm8aYybX6XDt5UXgOh1Yaioivk1CeDjB2FZ4nzwp4O3Nqs497l8yxg4Zh4QaV7QSY+wlxlgmYywzNzfX7LW8HhL2PK6H718SUfqSHj7DJYwfD/ToAQFvoeqo2t7M6Sphdw89BiuUnHY+UhLVpYScoep5U5KQN1mPFkmJODNCDz5UQlAQxbaq4H3rFi0X8PZ+OQpuA4AOAHoCyAHw94pW5Jx/yDmP4ZzHNGvWzOy1sMMKBh6kBKDifxtQskWBjw+cBu9BgwS866Cc4s1WpxT03W/ALomS0yLOKTh1yqX7LVT7VX1vKgrC1hpw9vl4hG8wIP2vNFOYLfBetsw2eHfpQvDOyKB1Bbw9Tw6Bm3N+jXNeyjkvA/ARgL52b8RYvs9vvYyG/0rAmokySibpnArv4cMFvOuanOVNNkWH8+/I2PYUJadNTiFvCgk5qmp70xgzmSyj/coEHPuzjB5v6crBu2XL6sH7uecI3ps2lYd3YCDBOyfHWd+KkCNyCNyMsZaafycAOFLRuhVKU76vZ0/gidckfD9gIa79bgmKi+EwvPPy7IO3eM67dsmZ3oycI4Ex4HykhF2DF6LLBlHyVMhxVdubmpjJGBAzX8KlGQvR5L9LsGkTHsF7xgzXwXv2bIqfy5cLeLtTtjwOtgrAXgBdGGPZjLEXALzLGDvMGDsEQALwG7s/WVO+DwB63lQwPGMxtvaajzVrqESfrfC+dImWd+hARVrsgXerVgLe3ipXezNor4KwMCp5OnjXYuwaMB+lpc49BqHaKZd40yJmsu0KuqxfjNsvzse+fRDwrkNivAYfc4mJieGZmZmP/i/dqlBW+Yt6hCQbAFnGD40kbNhg6kH7+QFlZUBqKnD4MDBiBEEYAAoKyIB375JR27Sh5WfPAqtWAWFhZNSQEDL0tm3A7t1ATAxVE2KMTPzppzSt3aRJQFRUjX0ddUaMsSzOeYy796MyWXrz3lcK/KbrcPZpPdptMmDtczLOR0qIiwMiIty3n0LOldd5U1HAdTqclPToohjAZBl8qIQtW4C9e03FVhgDHj6k6pM5OQThrl1pE1eu0HB3YCBBuFEjWp6RQZDu0oXW9/UliKekENRHjybIAwT9Zcsofs6aRRcFQs5VZd50a+W0O30kZMXqEfL3RNyPo5KnvXoB48YRfFevNu95d+8ObN1K8AWoxz17NlCvHhnUWs97xYryPe/MTGDjRus9b7WQgVDdVtEgCem99IhaZ8oqB4Aj9g+8Cwk5T5KEW9P06JqSiP399HjQn4bNn34aGDAANvW8w8OBmTMr73mvXSt63p4st4K78QEF/Q8YsGd4PPCBAbdTaQioInhPmGA/vHNz7YP32rUC3kLkzX4/GLDjKVNWOQCcOePmHROq21IUNF5lQP4v4tFVMWDr6woePEC14G3tUbETJwS8PVnuA7cxQ9InRUbHzxKwYboM/xk6p8N76lTr8B44UMBbqAIZvVm0wlTydHKKDhHnFNy5Y1a8Skio5mT0JWQZTf+TgJsfyBj2ga4cvPv3J3h/841t8H74sDy8R40S8PZkuQ/cmgzJ5s2BYYkS0qSFyF+4BPn5tIq98I6LKw/vjh2tw3vECAFvoQpk9Gb9cRLat6es8r1PLcSA3UvAOXD9urt3UKhOShMzAaD1TAn3f70QXb9cghUr8AjeP/kJwTsjo2J4nzxJm6wI3v36VQzvzp0FvN0t94HbmCFZupV62M2PKhi6dzEynpyP5GQ4BO8GDWyDt2rwyuAdHi7gXWel8eYTT1BW+cAdi7F30HwANJGNkFCNS5NVXlICQFHQ7OPFCI6fj+vXYRe8ZdlxeOt0At7ulvvALUm4+18ZReN1uP7yokfD5sMSJZSWwqXwXr68anjPmCHgXWclSbjxgYzC8TqE/ZtqladMlh8lqIn73EJukSQBsozSyTpkjlmE0sk0bN56pgSdDm6D9759tK62SIuAt2vl1uS0gJESTg7Xo3lSIq5P0j8aNo+LQ6Xw1j7nLeAt5AqF/lTCoYF6tPzIPKscAG7cEPe5hdwkSULxC3r035KItJ56XO5MvuzcGW6D98aNJng3bizgXRNyL7h3K3hijwEHxsUjdIUB5/5rHDavAt5nzlQM7z17aF0Bb6HqKHCPgphMqqMfk2lAjzzzcqenT7tpx4TqthQFQf8z4OH8ePRON+D7N5RHk99o4b1yZfXg/eCBffe8BbxrVm7PKmeyjKiUBOz8pYzmv9JVG95btgh4C1VTRm/6rpVx6UXKKh+TrHv0SBhAHhESqlFpssqD3k1A6SoZE1bprML72jXb4f3YY+XhPWuWbfD286sa3qK2ufPlEVnlAQGUVX507EKU/HXJoyIXjsA7Otr18FYfpxCqpTJ6kw2TMGYMZZVnDKOsclXnzrlx/4TqpiyyyuuNlYCF5MsVK2ATvPv1Kw/vmTNdC++AAAFvZ8vtWeV5KTSzTcBuBbHbFuPshPlYvx42wXvs2PLwnjixanhnZ9PyyuA9YEDF8NY+CylUC6XxZuPGQPsLCvpsMWWVA+S3vDw37qNQ3ZMmq/zKFQCKguCli9HsnfkICUGl8H74kGLbyJEmeH/7rePwHjlSwNudcmtW+a0PZQTP1uHU1EXgxmHzYYkS2rZFpfC+cYOW9+7tGLxXrLAOb21Sh1qFSMC7DkqScO3f5M3LP1+ESbIOXzxPWeUhIabV0tLct4tCdVDGrPKSSTqcnrYIxRNp2LzeWKqhXxG8r16l5ZbwTk+vGN7q3PMV3fPu398E73Xr7IP31as1/9XVNrk1Oa3heAnXJ+rRRU7EqWF68KE0bP7883gE76NHad3mzenqr7SU6uu6At6WGZlaeKslBAW864aaT5FwbqQerf+biEOD9CgeTMOTPpozRvWmkFCNSZLg84oeT+1IxJ7H9UgLJl82bIgK4T1liv3wXrPGBO9Wrczhffs2LVfhffy47fD296fRTQHv6smt4GbbFUR8Y8C5GfFo/ZUBe942Dpsb4d2mDRlCDZAtWrgP3tr6v9oiLQLetVNsu4LonQak/SQePXYZ0PKEgpAQmolO1cOHpls3QkI1IkWBT5IBZX+KR78fDDj5gfJo5EeFd3Cwa+G9bJl98FYTORs3piItAt7Vl/vAPXcuMGECmCwjYnkCTibI6PPWBFwYNfcRvKdPdy28Q0PL3/OeMsU2eAcFCXjXWmm82foTyip/cukEjN80F4yZr/r99+7ZRaE6KKMvIcvweSsB/utlTEuZAL9fzDWD9+zZngPvTp2Ar78W8Ha23NrjVqtYMEbJZr5+HHn5wFdfodrwlmXr8N67l9Zt0IAMFBJiDu9OneyHt2XxfqFaIKM3W7emYMg5R1mZaajc15d+q4k8QkI1Ik3lH19fwN+Po3Fjgq6r4K29520vvHU6AW9XyH3gTkoC/zwVxRN1KPrDIrApOvhtSMWdd5Owf7998LaWsHb6tHV4b97sfHhbm3lHyIuVlISy9akonqBD6Z8WYfT/dFg/IxVrhiWhtJTavLSUVi0qAi5ccO/uCtURJSUBqangOh0uxBkTelNTEfFtErp2tR3e9iastWjhGniLe96Oy6097vzHJaT11CPgnUQU/1wPNkyCJAFDhsAueJeU2A7vqCgBb6GqdbmzhL099fD9SyLOj9bjShcJQUH0WkSE+brbttX47gnVVUkS8nV6tFueiKND9CgZIsHXF5g8GVXC+8oVWt6li33wnjGjcnjbkrBmDd5Nmgh4Oyq3gjvssIKBhwzYOTQeJf8x4OEmBYzBafD+6U/Lw3vSJOfAW53zVsC7dqrNGQUDDlLJ07YbKTlt6FB6zfL57UuXqOctJORyKQrCZAMuvxCPyG+p5GlJCSqFt5qwtny5bfDu25fgvXkzxbbg4Mrhff9+eXj/5CcC3q6U20ue+q6V0eqTBKzTyeA6nRm8Bw+2D97Fxebw7tPHPnirj1PYAm+1hKCAdy2U0Zv+62X4vp0AeZKMCat16HyZSp6q4FbvcwMU6ISEXCpNydNWHyfg0t9k9P+nrkJ4q55UZ+2qDN7aIi2jRhG809Ich/eAASZ4r18v4O1sVQluxth/GWPXGWNHNMuaMMa2MMZOG383tvuTNeX72rcHBv5Jwm5f9KEAACAASURBVK4hC3Ft/hLcv08GGjbMPnjHxVUP3tpnIe2Fd1UT1gs5XzXhzUGDgLKnJOwaTOV4AaBpU1qtUSPTW9QnFoSEABd506LkaVe9hPwXF6LNmiVISUE5eH/zjXV4Wxs2z8mxHd7Nm9sH72PHBLydLVt63MsAjLJY9gcA2zjnnQBsM/5vn4zl+/h31Itpf0HBsPTF2NlvPpYvR5Xw1j7nfewYbdJeeFu75+0ovG2ZeUfI6VoGF3uTMWBgoYLBuxbj627z4edHvvPxMQUogAKemKdbSKNlcLY3NSVPOQegKGj76WKU/XY+Tp2ijoI6a1dl8A4Kqh68Z850LryzsmhdAW/bVSW4Oec7ANywWPwsgGTj38kAxtv9yZKEwhUyCsfrcPvVRY+GzQf8UUJ+PqqEd2CgCd5r1zoGb19f2+CtZmR26lT9OW+FnCdXevNBsoyHz+pQ8JtF6PBHHdY+J6P0SQklJVQet1070wWgqi1bHDsOodonl3jTWPK07DkdDjy7CGXP0bB5V72E0aMJorJc8/BWp7h1FN5ffSXgba8cvcfdgnOulou/CqBFRSsyxl5ijGUyxjJzc3PNXisaJOHQID0avpeI28/rAUlChw5UwawyeH/9dc3C29rjFNWZ81bIpXKKN+/1lXBwgB71lybi8lg9zkdKGDaMPHf5smmYXPts9/XropKaUKWqvjclCfdm6tHry0Rk9NHjZk8aNu/bF26D95o15vCeMcO58L52zUnffi1StZPTOOccAK/k9Q855zGc85hmzZqZvVY/U0HMPgMyRsXD/2MDrq2mYfOq4J2VVR7erVuXh7etCWuW8LbMyKxq5h0Bb89UdbwZdlhB3/0G7Bsdj7C1BkSco+HJTp2ordUJcBo1Inir+vJLVxyJUG2Tw95UFNRfaUDBa/F4fLcB215XcPMmvaSFd3WGzZ97rnrwbt3aufBOThbwtpSj4L7GGGsJAMbf1+3egjFD0idFRjc5Ad++IKP+CzqH4T19enl4P/ZYxfB+5pmK4W3tcQoBb6+RU73ZfV0Ctr4oY3KKDldX0TSfnFNAAYD27em32uu+cAEoKHDGYQjVQlXPm5qs8vr/TMCDZTLGJOuswvvkScfh3bVrxfCOjaXYuGWLgLc75Si4NwCIM/4dB+ALu7egyZCsXx8Y8baErJELUfDnJbh4kVZxJbxjYpwPb9XgFU1YL+BdI3KqN4ODgX5/oKzypv9b8ih4jDfenTx2jKAtet1CNqh63rTIKm86WULx7xai17YlZlNuuhLeo0cTvPfurR68n37aOrw7dhTwtkW2PA62CsBeAF0YY9mMsRcA/BXA04yx0wBGGP+3T8YMyQcbqYddP1PBoB2LcXjkfHz6KVwKb/XqVAtv7eMUjsLbWhUia/C2nLBeyDHVlDcbZFFW+eln5z/Kog0MJE8UFlJQ0c7Tffo0cOdO9Y5NyLvlEm9qssrv3gWgKGj4/mI0ems+Cgup/Kg98M7IoHVrEt7qeTFwoDm8y8oI3lOmVAxvPz9xz1sV45qi9a5WTEwMz1Qf3ANwO1WB/wwd8ibr0fZrAyDLKIiRkJxMw43Tp9O83ABw9iywejU9QztrFgVKzqnc5O7dpuFvxiiYfvopPco1eTJBGKAMxeXLKdDOnk31cgF6lvDrrwnCzz1HBiktJUMdO0YQ7t/fuM+3yYD375NRW7Wi5WpSSIsWtDwoiPZPLYTQrx9thzE6AdQJ5XU6uq9Um8UYy+Kcx7h7PyqTpTfz1yoIijN6c6MBy5+R0eHnEnJzgUOHKCDVr08XmPfu0XsCAkwV1CIjyadCni2v86aioHSSDnue0GPAQQP81smAJCEnh2JbYCDFNjV5MiODCkSp9659fSm2qTUmRo8myAME/WXLKH7OmkWjgwAVkkpJof9nzDDFtk2bqIaF2oNmjG4ZrlhBiZpTptDwN0CxeOVKituzZ1PBK4DqH2zZQjF60iQavSopIfifOUP5SH360Lo3btD+lZbS/rWoMLWvdqgyb7q15Gn9cRIujNajbXIiLj6jfzRsHhdHQdFazzsvz7znPXw4MGhQxT1v7XPe2p73smXle95qUoe2592tm3097+rMvCPkOao/TsLFMVQT+vAgyiovLqYABVAgyskhaKuBLzjY9P5z5+ixMSEhp0qSwF/WY8j2ROx5XI8TLWnYvGVLim229ryfe45gvmmTec979myKn8uXk78B5/e8ly2z3vNet67qnvfs2bT/db3n7VZw+3yvoOt2A45PjkfTFAMOv2ccNreA96VLtH6HDsC0afbBu1Wr6sF70qSK4e2KOW+FPEMBu8mbOT+PR4ctlFWem0tzuPv4UE9CvWcXGUm/tQVZAPKRkJBTpSjw+8iAkoXxiM0yIOMd5VFpZW+B97179sN7/35aV8Cb5D5wGyeFZ7KMLmsSkPV7GZ1+PwFXn50LwATvevXIKJ4Ib1dNWC/kZmm82fKjBBQulzFl9QS0f2cudu4kTwYHA2PG0OpffQWEhZFnAVPGeW6uaFMhJ8roS8gy/P6SAN+1MqasmYB7M+aawXvmTNfD+9NP7YO3WlWwdWtaXhm8rd3z/vJLAW+t3NrjVieF9/GhRDMfX47sy3TPGqBAOHt25fBescI6vDdutA7v48dpG1UlrFUGbzUjU8C7FkuT+9GkCcAYR0AgPZFTWEgBsVcvU07F7dsUjJo2JU+p+vxz84xzIaFqSePLgAAgwJ8jJMR8UqPw8KrhvXZt9eB95Yp98F692jq8rSWsHT1qG7zj4uouvN0HbuOk8CWTdCh7fRF8purgtyEVFxYmYetW2+Gdm2sd3pmZ1uG9dm15eBcVlYf3mDEVw1v7OIWz4S3ueXuAkpLAPzd5EzodNr+SisO/SMLYseSXS5foJyyMfOXnR4GmY0fzTT18SLAXEqq2jDGT63S4/vIicJ0OLDUVkZuTys1IqIW3tUfFTpxwHrwLC03wjokheG/dag7vZs2sw/vuXcfh3bRp3YW3W3vceT0k7O6hh8/biSibq4fPcAkTJgDdu8NmeE+d6hp4x8Y6Dm/1cQpHEtYsp80Tco8udjB58/bzeuT3kFBURPO8d+9ObbZ8ObXn7dt07w6gLNuAAPKbql27xONhQk6SJCFvsh7NkxJx9mk9+FCpwumEVXg/fFge3qNGOQ/eK1ea4D1mDMF7zx5zeM+aZR+8R4wQ8K5MbgV32GEFAw8Z8P2T8Sh6z4DSrQp8fFAhvK3d8+7Y0TnwnjmzanirBq8I3vbMeaudsF7A2/PU7kcFgw4bkPZ0PPw+NiA0Q8HDh/Sa+hhKt2702MutW/SUQaNGpkfCiovN5+v+9NOaPwahWihFQdhaA84+H4+WX1ByWmXTCVcE7379nAPvyZOrD+8ZM8rDe9CgiuHdoUPl8L5uf61Er5P7wG0s3+e/XkaDpQlYM1FG8URdpfBu0MA2eKvlR4cPp6s3a/AODzeHt5rUURm8tcX7rcG7qgnrK6r/m55uqv+rnoAC3m6U0Zt+62T0/ioBx9+g0pKhGQoOHDAloT31lOkZ048+omFzgLxUVmaatxugYKIGQiEhh2T0JZNltF+ZgGN/ltE9UWcV3toCT66Ed7du1Yd3mzb2wXvq1MrhnZxc++HtPnBryvf16gU8/qqE7wcsxNV5Sx5Nl+govJcvN8F7xAjr8J4xw7PgbW3CegFvN0njzYAAIGa+hFOTFmLA7iX44gsKRgAFmn796G/GKBA9fEheYYyCB2OmzX7zDb1HSMghaXzJGPny0oyFaPzJkipLKzsC786dTUVWAFNsE/B2v9wHbk35PgDodUvBsIzF2NprPtasgdvhrSZ12ALviur/Wt7zFvD2Ell4E4qC6A2LkT5kPsaONQW9tDTqffv4kEfVKlEbNpAX/f2pLVV4c06eqsFihUK1SRa+ZNsVdFm/GLdfnG/TvAj2wlunI3hv3GiCd+PG7oW3+pSGJbx/+IHWrSvwdh+4JQmlq2QUjtfh/u8WPRo2f/xVCWfOoEJ4q70dV8NbLWRgyz1vV855O2OGacJ6Ae8aksabBb8hbx77s4yzbSX07Am8/DKtdvo0JfrWq0cemTqV7nHn5ZEfi4vJF9rZxPLyRJa5kIOSJECWwXU6nNAZs8plGTHzJZtnJFTh/eCBc+C9YkXNwvvIEevw3rChbsHbrclpt3tLyIrVI+Tvibgfp380bD5uHCqE95Yt7oG3avDY2OrPvOPohPWybKpCJORaXYuSkBmjR/2liTjypB5Xu1FpyaIiategIGpHPz8KLufOUbBp395UoAWgZMl69cyf7d650/TYoJCQXZIk3JqmR9eURPzQT48H/aUqpxO2nNQoPJximzPgHRAg4O0OuRXcTQ4q6H/AgD3D44EPDLidahw2rwTe0dGuh7ejJQTtmXnHUXhrSwgKuU7hJ+mJhwuz4tH+W8OjeeKzs+l1dZKEl1+mkqfFxcD771N7FRQAzz9P5VEvXqS2B8zvdycnk6+EhOySoqDxKgPyX4lHF8WAba8rj2JbRdMJa2tEOApv9Z63s+C9bZuAd3Xk9qxynxQZHT9LwIbpMvxn6MzgPXZseXhPnFg1vNXg6ii8tc9COgpvy4zMiuBd0YT1At5ulCZ7t11yAvzWyZiyXoeIcwo+/ZS+/4AAArSvLwUjgPIR1MB47BgwZAj9rd5qeewx00cUFwOffCLudwvZIaMvIcto+n4CbhhkSB/osO11xep0whUVeLIF3iNHmsPbz882eF+9Ssu18NYWaVHhvXu3Oby1RVrOnqVtVAbv4cPtg7ePT+2Ct0dklTdvDgxLlJAmLUT+wiW4cYNW6d3bMXivWGEd3tpHxTwF3hUV7x81ioblBbzdII03ASBgpITbr1BWee/eNCyenU3FHq5epYADkF8nT6a/FcU0HN6nDwW2nBzTI2MAeTI1tQaPS8i7ZeHLNrMk3Pv1QnT5conD1Rkt73mrE+X072+C97p1tsHb3586SFp4T5pE54ElvPv0MYd3SIgJ3qtWVQzvggJaPniwdXirz3lbwnv27NoFb7dnlZdupR5286MKhu5ZjIwn52PZMrgE3tevl4f3gAEE702byj/n7W54q/V/09KqnnlHyIkyepN/Z8oqD/toMfYOmo/u3YHXXgPataO2Tkqi4MMYJZ5FR9NrISEUVHx8qI1UoOfl0W8f45l36JApJ0JIqFJpsspLSgAoCpp/vBhBr8+v1rwIrVqZ4L1sWXl4Hz9uG7xnzy4Pb3WebUt4P/OMY/BetqxyePv71w14uzWr/O5/ZRSN1+G6njJ3fdbKkBIklJaixuD99NME7337TPAOCnI9vC0zMp0xbZ6QkyRJeJAs48E4HX6csQhlz+lw8wMZ5yMlFBaSP7p3p1X79wfOn6d2OXiQAla7duQv1Yc3b1JCWqNGVGGPMfOJR775BvjxR7ccqZA3yZhVXjpZh8wxi1A6mYbN28ySqj2pUXXgnZlJ6wp415zcmpwWMFLCyWF6NP8gEdcnU1Z5ixYER2fAOzS0/D3vKVM8A97OnjZPwNu5etBfwrmRerT/NBE7o/X44g4NT6plT9XqadHR1ANv1oyCSlKSCeQlJZSkBlCb3r5NP8OHl/+8lSvrTp1loWpIklD8gh79tyQiracelzuTL50xI6Gj8P76a9vhnZ1tO7zDwsrf854+vXJ4p6aaw7t9+4rh7c3lUd0L7t0KnthrwIGx8QhdbsC5/9HQpC3wlmXr8N67l9Zt0IAaKCTEHN6dOtkPb8vi/TUBb3unzRPwdq6aHFQQvcOAwt/HY8BBAxpkkTc3baLHWVQVFJBPoqPp/8GDTQBWM3sbN6aazN260fKdO4GICPMsc86pbKrqJSEhq1IUBP3PgIfz49Er3YDv31CcOiOhq+E9ebLt8J41i0CrhXfbtiZ4W7vnffiwObynTrUO77g4+mxvhXe1wM0YO88YO8wYO8AYy7TrzZrM3ai1CdjxCxnNf6mzCu/k5PLwPn3aOrw3b3Y+vK3NvFPRc94VzXkr4F2zcoY3IcsIfCcBAZ/LmJxCWeUNG9LIzurVtOqpU9TGaoJaVBTw6qs0RJ6fTz1wxoALF+gRx/btyTNqr1ybaV5aChgMoixqbZfD3tT4MujdBJStkjF+la4cvKs7I6EW3tYS1qzBu1On8vC2lrDmTHgXFNgP7wMHaN2wMO+GtzN63BLnvCfnPMaud1nUgx7+loQjYxeiZPESHD1Kq6jwLikpD++f/tQ6vKOiah7e1qbNszbn7eTJzpnztip4q89CClXfm7QVCfgDZZV360ZD40OH0ksHDgD//KcpqFy/Tu0RHU2eHDKEAgvn5LMePWi9Ll2oPa9eJR+pKioC/v1v4N696h24kMfLfm9a+LLeWAlsIfnSWo0Ia/C2NiOhtdLKKrzv37cN3jpdeXg3aeI4vL/7zjXw/uKLiuGdm+toU9a83J5Vnr+WZrYJ2K2g77bFODN+PtatQ5Xw7tPHOrwnTfJseDtjztvK4G1ZyEDIAc2fD754MW6sI29CUcD+uhj7npqPwkLy01NP0e+ICBoGV4fhtm+nZ7hbt6Y279yZeuB+ftSuX3xBj4bdvGmawzsnh7yhTgNaVAS8954pGAkJATDLKs/JAaAoCF66GM3emV/pdMJaeFubkbCieREqg/dPfuJaeO/aVT14DxtmP7yTk70H3tUFNwewmTGWxRh7ya53ShJufSgjKE6H09NMdXeHvyWhTRtYhXdxcfXgHRfnOnhXNGG9K+EdE1Me3taqENVRVcubNz+QEThLh/RRi1A0Xoecf8nI6SqhqMi0WoMG1BZTp1IvPCSEhrlTUuj+NkC3TEJDqaft50dBpbSU/LZ7N4Hdx4cCS2kptT9g6nmLe961Uo5505hVXjJJh1NTF6F4Ig2b1xsrVToj4dWrjk9qVBG8BwwwwXv9eufDu3fv6sF7yJDaDe/qgnsw57w3gNEAfsEYe9JyBcbYS4yxTMZYZq7FN9JwvITrE/XovCYRp4frwYfSsPn06bAK77i46sG7YUP74d2/P8HbWvF+Wyas79u36jlvtdPm2QNvtQrR3r1V1/+tg6qWN0OekXBnuh79NycivbceH56WcP8+FV85cYJ82KCBqZpTgwZU+rRePWDaNAp4AAWflSsJ3kVF1EN/5RVqv/PnKUCWlVEwadqU2lZ9xru4GPjPf0S2eS2U496UJPi8osdTOxKx53E90kNo2NyV0wlXBe9jx5wP75/+VMC7MlUL3Jzzy8bf1wF8DqCvlXU+5JzHcM5jmqkZPEax7QoivjHg3Ix4tPrSgL1/MQ6bB9BjNJ4A78qK9ztrwnptIQNH4W1L8f66pOp6M2ivgpafG4D4eAw+YsDsdgqCggjUa9YA775LbXPrFgU0gL7vW7cI4NOmkff8/CgI7NpF62zbRt7s0YM8MHAg/X/sGHnTz49+VKlFXkTSYe1RtbypKPBJMqDsT/Ho94MBJwyKS6YTtpzUqDrwzsqidV0Bb7X+gQrvO3fqBrwdBjdjLJQxVl/9G8BPAByxeQNz5wITJoDJMiKWJ+BkgozeiRNwYdTcRxXMtPA+dozeZi+8rSWsWcJbzcjs1IkMd/06LXcnvK3Nedu1a/Vn3qkLcpY3IctAQgLVLH9tAp7dOPdRAOvVi4JYaSnwt7+Zt7VaHa1DB/Lh9OkUKOrVo7b8178o47y4mLJv58yh9c+fp/WLisj3qjgHPvvMVKdAyHtVLW9qfOnzVgL818uYljIBvq/MddmMhNonVRyF91dfVQxvdTRJC+/PPrMd3qtWmcN7xoy6Ae/q9LhbANjFGDsIIAPA15zzb+zagnGGBcYoEPr6ceTl01WaJbzXrnUM3r6+VcN7xQqUe5zi2jX3wtuVE9bXATnNm9r/fXzJc+3b0/c9Zgy91Ls3ZYGrJ/+aNcD339PwOEDBqEsXuvUCkNfVIfZNm2hIMiqK2jE2ln5fukRBTqstW8g3YmISr1b1vKlpfF9fwN+Po3Fj104nbAnvGTOcB+/kZHN4T5pE3hfwrlwOg5tz/iPn/AnjTzTn/G27NpCUBP55Koon6lC8cBHYFB38NqTi9jtJyMoqD+/WrcvD29aENUt4p6XRugLetVPO8GbpOvLmndcocRKpqTj0SlK55DSAfPXKK6Z712VllF2+ejX9n5ZG03t26ED/h4ZSMlvv3rTujh3k69JS8olabU2bmKYWazl6lO57P3jg+Pcj5D5Vy5tJSUBqKrhOh4uzjQm9qamI+Dbp0XTCtsC7ujMStm7tOnhHRwt42yLfN954o8Y+7MMPP3zjpZdMSZT59SNxYOddRK5MRPGv58H3xTmIjKQGT0+nDN1Onajxo6KoiEVaGkGoWTMaeuzQgR7FOXyYTBkcTNCsV4/WvXqVQOfnR0PNeXm0PCiIDBgURMuPHSNTRUZSQG7alIpjZGSQIaKjqZE7dCAzpadT8OzYkZZHRVHiUno6vS8sjMpitm8P7N9PAbdrV9PnBgfTfly/Tvvn60u/r12j5aGhdHUbHEzHdeQIHWf79rTdZs3ohEpPp8+NjqZj7NSJen/p6abeob8/vX7mDB1PeHj53pwr9eabb+a88cYbH9bcJ9ovS29eCYjEod130UVOxK6+87A9cg7u3KFgEBtL33VZGX2fkZHU5qGh1M7NmlGeQtOmdDF46xYFhAMH6CLs+nUKkhER1E7R0eSNK1coQJw7R/65c4deU4OGjw8FrAcPyCPt2lGAFnJcXufNyEjkX7iLlh8n4uioeWj6uzmPYsf16+SL4GDz2Hb0KMUgNXaEhVHHR40dUVEUIzp2JBhnZNDtGjV2REVRIlhGBnVQmjY1PQqZlUWxU41tbdqQx9PSKNZ27UrnSrduFIvT0mgfwsNNse3wYdq/jh0pbjdvTvuoXvCquSKdOxOE09Pp3IuIoHyoqCi6qMjIoONu3Jg6Ze3aUSfoxAn6/MBAWubrS9u4eZM+X+VLdjZ9ZuPGdD6HhNBnHjxI526nTqZRtJpQZd50K7hD0hW0fm8Bdvefh8dSDSjtFQv/TpEeCe+zZ22D948/Vh/eqsHtgff58yaDexq8vS44AmiQpaDd+wuQFzcPHTYbcKVVLE4XR6K0lB7jOnmSehpXrlA7t29P77twgS6+hgyh75gx8s4zzxB4r14lz+zdS/e5AwKo7aZOpYlL9u2jgKHC+s4d6slcvGi+v5xTMFH9py2fKmS7vM6bioKQNxbg8tR5aJlqwO7CWLR9KrIcvENCKHY4Cu/0dPvhrcKxTRvydXo6xVprsc0S3ocOkZ8rgrfaMencmc4Ja/A+dYrOHy2827alzPaTJ83h7eND27h1y3Ph7ZngNpbv81kro2TmHGy+EYtub+icAu+uXcvD+9q1moF3dHTF8M7KohOoWzfzz01PN4d3VJTr4B0VVfPw9sbgCJ0OLEVG6C/nwLd/LLq8rkN9KRaniiMxZAh9t+qjXJcuke+uXqVgl51N2eK+vhQg9u8HnniCirZ0724Kftevkw/LyqhXEBpKbZSTQwlr9++Tby9epIBTUkKe0pZEvXyZgop6QShkn7zKm5qSpw1enYOzjWMRnaDzGHirsc0V8L5woWp4R0eXh3ejRgTvffuqD+8DB+hcqyl4eya4f/ELuonw/PNo0gQI7R6Jg0f9ELR2JQJmT0dAAMGzpMSz4f3wIS1/+JD+rw68c3Pp8y3hXa+eucFVeHfoQK/ZAu+SEjoe9eq0JuHtVcERMPMmAPri/PwQtG4lMjpOx7RplEswaBB5TQ2CP/5oesQlK4ugyjl5NSSE2iM4mIKLjw/dE+/QgaB9/z4FltxcCka5uVRaMi+PvN+qFQWYu3fps9RH0ADyXno6bbNtW9H7tkde5U0LXzaLjUROrh9CP1+JbS2mW73lZgnvI0c8A945OQLeVckzwR0eDrz6KnhMLFhkJJocVBDxj1fx1fClyMyLRFQUagTeubnW4X30qHV4p6eXv+etBs6K4K0avCp4p6XZDu/Dh22H9927lcO7VSvXwturgiPwyJuIjaUvTFGAV1/Fxd8uxZG7kejbl058Hx+6t1ZWRvF04EBq6yNH6Hd+Pl0Eck6B6uJFChB+fhRc+vShdvP1pe3odBREr12j96rFKx48oPaZOJF8o+1x+/ub5vY+f94UVEJCavTr81p5lTc1vuQRkWDbFTRa9CquL1yKHRcjH8U2W+HdoYPz4Z2ZaduweVSUCd4NGtB2ahLeKl88Gd6eCe7ISBQ9EYuSSTrcu3YXQYsWwGctle/LyKBA5k54d+tmHd4tW9YsvLUnoBbeXbtah3ezZvbDOz3dtfD2quAIAJGReNA9FqWTdTh/+C7qv70A+f8nI7+H9Ai4KhjPnaMANHAg9XQbN6bM1+7dqWPUuzcFhLw88sTRo6Ze+fHjFHgaNya/BwUBI0ZQm2dkUCAKCaH337xJ67drR1Bv3Zra0HImMXUEqKCAfKFWYROyLq/yZmQkEBuLsud0OLT7Lpr/YwFYioxmOgkhIeaxrTrwVmOHrfBWR+3U6oGVwTs/37xj4ip4W7vnrYX3iRPOgXfnzq67SPZMcAO43yISR9PvosOnibj94jwE6eegSRMymL3wTk/3LnhHRtoGbzUj05XwVjMyXQVvrwqORt1qFIkzB+6i+3rKKl9bf86jEpA5OdQ+BQWUoHbxIvDkkwRuHx9q07Iyqo4WGEhtdegQ1c4YPZqCx5kz1A6XLpl65Veu0PZ8fGhb584BP/sZJacdOkQ5Hzdv0rbv3KFA1KwZLVOz3FVphyK1U4cKmcvrvBkZibs5d9FueSLSh8xDvV/NQVAQnbvugnd6unV4Hz9eHt5pabbBu3Nn6/Bu2pQ+z93w7tTJlLDmKnh7LLgD9yhouXQB9j05Dy0/N+Bmx1iEdo+0Cd737pnD+/z5iuF95Ih5tnloqPvhrRq8LsDb64IjgNAMBeFLF4D/dh7abjSgw5RY+LSPUtyFJgAAGudJREFURE4OfZfnzlFgUnvPP/xgusddUEDBqU8fCk6hoZSJ3qgRtUHTphQIL10CfvMb8qavLwUwdXg9L49gfuqUyfsnTtB99YEDKXDl5xO01WfHH3uMhtXVGh1lZTQseOAA+aVhQzd8sR4ur/OmoiAwfgEKXqKs8m/zYxE+KLIcvCt6zLQ68L53r+JHxZwNb7Vn60x4t2njPHiHhroe3p4Jbk3mboPXKKs86g2dVXifOVM9eO/fbw7vVq2qhrf2WUgtvNWkDkfhnZFBDW8Jb+2zkM6Etzq0ZCu89+1zPry9MTiq2btszhyw2Fg0fFGHpiNjsS83EuPGAePHU6Z4QAD1uNu0oe/19GkKIEVFNGR+4AC9/vAhBYPWrSmIBQWRL8PDqU06dSIP1K9PJVKbN6e2v3mT/H/6NAWVCxeoffv2pYuH4GDg8cdp3Tt3CNqWve/CQtqPkycpONXks6ieLq/ypsaXgS/Pwe1OsYh6U+dyeFvWiMjIqPpJlZqCt+Vz3iq8OS8P74wM++B9+zZt0xLeTZrQ96OFtyuGzT0T3JoMycBAoPWQSOw/5IcAeSXujZ+Ohg3pCwoPdwzezZs7F97ahDUtvNWENW1SR2XwVu8LORveFSWseQK8vSo4AhVmlfutXold7aajc2dqv+Bg+snKonKmo0fTtJ0hIeTXmBjy2I0bFATu3ycf7tpFPi0qorYMDCTQBgZSAOjTh9qzc2fySq9eNFzu70/BODubvFxaSv4qKqLPzs+nz2jYkJYD5hC/e5fa9vRpClAigc3LvGnhy5DoSNwv9EPDL1cixX/6o9ihPj7qLHg7WuBJm7BmCW9/f9vueVvCW00IU+FtrUiLM+DNmHV4X7pkHd6uuOftmeA2Zkg+6E7PbQfuUdDmb6/i+/FLoZyLREQEHIa3WsFMC+/27csPmzsD3mFh9sFbLdJiK7wDA80fFfNGeHtVcASA8HDwV18F7xML1t6UVV74zlLsuRKJDh3oOwfoBN+7l4JB69b0f2AgfX99+9IUyrGx1M7HjwPDh5uqN6mV2E6coKCplmw8c4baprSUfk6cIDA/8QS185kzBPI2bSiw3LxJvrl3j9r24UP67KKi8slrAA3l79tH+9OiRd0eQvcqb2qyyu81j0TAbgVBv38VZX9fin25kXZXZ7QF3tWtztiwoXV4t23rGLx/+ME2eN++XT14R0RYh3d0dM3B2zPBHRmJ251j4TNNh6un76Lh25Qh2WqGhOPHqaFtgXdxcXl4R0eXh7eazV0VvKOiyEDOhndhoWPwVksIqo9TWIO3tRKCngRvrwqOABAZiTudKXt3//d30XTJAuz6lYyT4RJycug7Cg6mnmxwMNWGb97cVD0tKIh61U2amJaFhBDgO3akqmpRUfSeI0eAsWMJyuq0oHfuUK/6+HHqrXNOATU7m7xZUEDPiE+cSElxFy7Qe2Ji6PVbt+j1+/dpX0tLKUAxRu2sSp0YZf9+8thjj9W9Z8C9ypvGrPLSSTrsU+4i/D16EifkGcnh0spaeKs1IrTwdkZpZXfAu0sX74e3Z4IbQEDnSJw9cBcdPkvEJd08NHxtDgIDqWGrgrd6T7l9e+fC++pV2+FtWcjA1fC2t/5vRfDWGrxTJwJBRgYFeFfA26uCo1FF4ZHIO3cXPVITcXz0POx/Yg6ys+k7ys2lkzMjgwDNOQXG8+dNj4fl51OiWHg43TMLDaU2KS2lR8UAam+1Ot6gQeT18HBqt2eeoUfD1Me/7tyhC4UzZ2i7JSV0IXD4ML3/zh3qeY8dS73to0epfdu2pYuze/foPerQuVr3HKCe+alTdCx37tA+BATUeBO4RV7nzchIlN2+i4gVidjVbx4ePj+nWqWV09OtV2d0B7xv3LAO74YNPRfe6en0flfA22PBzbYrCFuyACfGUK3yUw1i0bxfZI3Au7J73pbwtqd4f0XwfvDAOrxtveftSnhri/e7At5eFxxBTzw0XrwAmDcPLdYbEPtKLAbPjMSePQTekSPpe2rVivzBGH2XV6+SN4uK6DvNyqIe+c6d5NG8PNN96rw88uvFi9SO/v4UQI4epYuDIUNMPfmMDOqVx8WRr2/cIFC3aEEB5d49+swffiBPhoZSj/vBA7r/3qABPW4WHEx+Kiyk49T2sNVCMXv30rkRGEjbr829cK/zpqLAZ+EClLxKMXPzjVgEdYt0CN5VlVauaXinpZWH95UrNQvvNm1o/2yBt3rP21Xw9kxwz50LxMeDrVuHpr+bg/TSWPRYNAF5+86h3rSxLod3RQlrISG0rhbe9s68Yw3eHTvWXXh7XXA0ehPr1lHR8NhYegj73DlktBiLsDCa3vWxxyjoXL5MbarX073nIUPIi9nZNJzdrh15r7jYlHF+8SJ9rwUF1BPOyCC479tHPWIV8Hl55JuiIjoXuncnj3XpQoHB15d2d/BgGtm5coW836CB6R73qVO0nDE8mpa0e3fa54cPTcPpWj14QMFrxw4aiq9f3zTcXpvkVd7U+NLn53PA+8QiOn4CsneeQ8FTYwW8axm8K/Ome+sqGcfqfHwo8Pj4cGRfpmdeAQo+cXEEnJUr6QsCCIJTp1KvZPly03zZI0bQM66ZmTSvtTqf94wZZJy1a8ksABlg1iwKXsnJ1HsBaJhx9Ghq0JQU03zZkyfDrjlvHZmwXpbNJ6yfOdN5c94uX2565jgqio7H3jlvV682JVDVeqnjyBb/BwTAbE5ugIKaOucvQH5u3Zr+DgujOdhHjKBHyAC6L71wIbBgAV0XAOTpoUMpkDVvTstOn6Z5vb/8kgJlaSnw/vvAX/4CfPQR+U71WFaWaaj9xx/pM37zG7reAMjvAwdSUCoqooCizvddVET73Lgxed3ysM+dI3+//Tb56McfzR83E6pBaXwZEAAE+HOEhFDsOHmSloeHU+x4+JBih9rO/fpRDDpxgmKhGtt0OoLLxo104QiQF+LiKH4uX04+A8ifkydTrFu50hQ7xoyhHIs9e4CtW02xbdYsAv7q1QRxgDw4YwZdVCYn08UsQLeLRoygC4bPPyeP+flRrO/QAdiwgUANUDyKi6P9T06mixGALkgnTiRwf/YZeZsxYNw4oGdPuhDdvp3WDQ2l/WvShNY9d46WR0RQ4v7t27RtNcHzqafoHD14kPalrIzaYNo0ujhPTaWLCYCOOS6OvofkZLoAd6bc1+MeOxbo3x8lk3RAwV34LFwAtn4d9sb8Cunp9IWoPUxtz1stguKKnrf2cYqKqhDZM/OOoz1va0VanNXzVp+FtOx5awsZVNbz1s55a6u8qlcDAGPHorhPf0CnQ0HOXQS8vgD3lq9D6Su/wpEj9D316GF677Vr9N0MGECvqVKLPqiVy4KDaZnaVn5+1Iu4eJGC13PP0ff/+ON0IZmfD/z2t9Thj4qiZQUFFIDU5Lj796mHcuYMeaSggIJFVhb95OWRV69epZ59TAx54NIlgvWQIeSzvDzaluX1iq+vaRnnBIFDh2h04PhxCozNm5sftzfJq7xpjJlcp0PuubsIeXMB2Lp1qP/HX1U5qZE9CWvOmJFQTXa17Hmro3aWPe+oKMd73pZTbtrS8wbMe94nT9K52bYtdci0PW91YhLLnvedO67teXvmUDmAvPqRyNp+F5ErE1H2m3nw+fkcdO1KASstrXbAu6Li/fbAOyLCOry1xfsdmbC+WbOK6/9WBG/LEoK2yKuCo1F59SKxf8dddJETsSN2HlYFzcHu3XT1feMGfQeZmXRSXrlC7ZmTQ+18/jy1yaVLpiv+mzfpvbm5tJ6aoObrS+scOkTL1BM7NJTaPCyMvn91buGMDPo9bhw93x0dTfvRpQswaRK1Xf36NEwfGkr+KSyknzt3KAidPk2BtbSU9vXBAzofOKf11Gz4+/fpf20im1b37lHve9cu8kp2NvmnQYPyPXdPldd5M5KSJpt/kIizz85D43lznDIjoeWTKpbwVi/8tfC2tTqjtkiLdl4ELbyPHTOHt5+f9YQ1a/BWbxtZwrtJk/JFWqqCt3p+aeGtnVVMhTdgO7w7dnSsSIvHgjskXUGr9xZgV995aPmFASw2Fr4dIl0Kb2tFWuoyvCsr3m8N3tZm3qlKXhccQclp7QwLcGvOPHTaZkCbCbFoNTgSN28SyKKi6Pvz9yefqc9d37hB97zVYbebN2lo8uhROnFv36b23ruXgPf99/S9l5XROkePkgevXKFtnj9P/s3Opt5uYSFtr0kT0zCljw9tu0sX8nnHjrTs+HG6BTVuHPWsb98mX8TG0rBps2amaUNVOBcVkS9zc+m4ABO0g4PJA+py7f3ukhI6Z48epeNKS6PvoLCQ/BsY6Po2dERe501FQcibC3D2WYqZB/yoapqj8K6swJMW3trqjI7Au6JJjVwJ7xYtvBvengluY/k+37Uybjw7B9/diUW3N3Qug/f9+/TlqXBU4d2ihcng3gTviqbN80R4e2Nw9J2mg0+KjOBfzIFP31g01evQenwsLvlF4sEDukfXrRv1ktWs71Gj6J7y4MF0P+z6dfLfiy/S8Hb37uSDs2cpl6B7d/r+WrQgrwPUDiUl5FcVolevUkA+e9ZU1vT4cQpS6rzfALXrzp2mYFxWRuudPUvnhY8P+ezHH8nDzZrRvt+6RZ/fsyflObRsSRcKJSXkRbX3XVBAPqioB66VOiHKmTN0kbJzJwX9Cxfo3qta9tXdyW5e5U21TLQso/G8OTjgF4voBJ1L4V0TMxLWZnhrn/O2F94uS05jjI1ijJ1kjJ1hjP3Brjfv20cZFZKE3r2BHr+WsGaijEOf7ENJCZCQQEkG0dHAli2U9AAA//iHecJadjYtX7nSlLC2YoUpYW3XLrr3uG8fsGkTBZy//pXqQbdsSQloJ07QNpKSTAlry5aZkjo2bqSEtZMnTQlriYmmhLVNm6ihAWDpUlNSx4oVpqSO1avpHmZODvDpp6aEtbQ0U8Lali20f++8Y0pYW7PGlLD28ccEjPv3af/UpI7Nm+mRHzVhrawMeOst84S1/ftp3ffeo/3z86Okk2vXaHlKCg21XrpE+6cO8WZm0pDsrl0UNzgH3n2XttGkCbBqlSlhrQavAauUs7wJgH7LMrBvHwICKFlMq/r16fe//qX9fNOjWkuXktciIqit/f3JDwMH0qZHjqQe8IMHFBh/9jNKIP7VryhQ9exJbbpgAfDaa6ZJSSZNoouDZ5+lXjRAwfjQIfqsiAha7+pVOi+SkkylUI8eJd9+/TUdGkCB5oMP6Pi+/pra+tgx8gDn1N7bt5P/69Wj7asXbYpiKgGrlaLQ77Iy+i5OnAB+9zvgP/+hc2jxYkq4UxP4L140LxIDWPdVRV6zd7k75LA3Nb5kDIj9vYSji2Tc3roP335rim0zZxK0ZRmPZrT78ENa/uABxQ412fWbb8h/x49TG5SWUuxQ8y2+/to0L/y//gXMnl0+2VWWyYvZ2RQ7CguBN9+kWgR9+lCy8bZtptgxc6YpYe3sWdrGJ5+YJ6ypyZ5bt1K1wSNHKPlLjW3WEtb+/W/aPx8f2j81YW3dOrqgvniR4pUa2/bvp3Pr++9NCWtLllBsa9yYEtbOn6fly5ZRwtqtW7Tte/do+fbtdJF+4ACdN5xT8qiasPb553ROA8D//R/xpayMjlG9WHfEm4xXdelc0RsZ8wVwCsDTALIB7AMwjXN+rKL3xMTE8EzVBVa0fz8dfMeO1Iic00GuX0+B5umnKfOQc4JWcjJ9gTNm0BUR5wS5NWsIejNn0lVNWRlBce9eCnDPPEPrPnxIwM/JIaN262Z6lnX5coLv7NnUiJwTnDdtosA5bRotKy2lDE21LGW/fqYknmXLyMSzZtFVK+e0XkoK/T99uinJaNMmOi8HDKATiXM6yVasIANOmUInEud0gqxc+f/tnV9sXMUVxr8TK5YlO1GSxqCK2LUlLFDAiGIURST8ewjCQSElqjAG6kSx6EsfWpSXoAo5POSBB6I+tC9VG9mJkhAF10qQQJELRYCEIpIqqtKgFhcSWkiTFgdiCJDYPn04O7276zvJrveu787195OuvHv2evacO9+dMzM7s2uxbdlivU5V69yMjlqP7/HHzTY5Ga1W37DBbiRV68EODpr/fX12o6vadR4etuv51FPR92i/+qrdIPfeazeSql37PXusrN5eu5Hi5CQiJ1T17lkJbRZUQ5uOo0dtZF0c54svAtu3F9pPn7a63rGj0L5nT7Ti1HH+vCXN4nNHRqzD+NxzkX183BKfWyXs7MPD1gA//3xkO3PG3q+z0zrCTleDgzYi3rTJdP/xx9bAjo1Z3ff3W9I9c8ZGyfX1ppGtW4GXXooa1ro609COHf4EW2wv5Vy3L76hAdi2ze6x5mYb/SxdaqOheK2Vaw9bm6qmyWPHTA/d3VHbtnevJdeeHmuzVG12Zu9eazs2b7brqWoDh6NHTQs9PVHbceiQJf9HHrG207UdQ0M2m9TXZ51S18l75RXrPPb3m03Vkv+JE9Z2r1tntsuXTZeff25J+Oabzf7JJ5b8Fy0y/xYvNvu771ry7+y0ToKqvf/Bg5b8H33UZrFUbeZxaMge9/VZB1rVEujIiCXU3l5r26amrG07edIS8IMPRm2bW5H/5JM2Y6Bq98O+fdGq+6Yms7/1lnUA7rzTdo+4j53277eYHnvMFp26L2saGrIOxpYtNitSrjYrGXGvAjCmqh+p6hUALwPYWEF5uOsuSy5u28DkpAW3aZMlo9HR6Fy3Vayx0ZKYo6PDhHfhggkUsJt23bpo5A3YhWposKTvRt6O/K1ig4N5Aa+KRt7AzK1ir78enbtkiVWK207huPVW6yR89pkJwPnX3W03xnvvRf7lbxU7eDAqY8WKwpG34557opE3EG2ncCPv/JHismXmX11doX+33VY48nb+bdhgI+933on8y99OceAAaonEtelwn9UWb4davHjmuc3N8WW0tNhf9yUogNVx3EcOXV2F5wF2ve+4IxoJOR56KFrd7RqCtjbbwuK2qQCmq6eftvpzmnBbYO67L9p2uWaNaWTrVrtX3JT8s8/aV2bfcovdAy7OZ56xhtX9nnh7e1S288tdv0WLCqfJi6fMXcPsOginTllHYmQE2L3bbC+8YDNzO3dax2nXLrMfOgQcPmyjyTffjLaX1giJalMkmrFx076ubXMj7/y2o3ibqWP16mjkDUTbTPNH3o7ibaYOt83UzYK6NRj5I2/nX/E2U0drqw0YJiYK/Vu71gYMbvQ6PW3v39MTjbwdy5ebfyKF/nV2WgI9ezZqrxYssLbNjbwdjY1Rx2b//sje1mb+XbxY6N8DD0Qjbxdjfb3dU62tplvHDTdY2dPThe13WajqrA4APwbwu7znPwHw65jzfgrgOIDjra2tei0GBlwfrfC4//5k7LVSxnzwb2AgqlcAx2ers6xps1bqp5plrFlTetmrVsXbb7+99DLa28uzU5u1rZ9Qy06ijFK1WXUB5h9dXV3XFGA+wEzb1FS8/euv4+1ffTXTPj0df+6338bbJybi7XG2ycnyyvjyy9L9u3w53v7FF6X7d/Vq6X6oqo6Pl+5f3LWOfKnNxjH/qFSbvusSpwmfjr/7Lt7+zTel695Xtk8/ly6Vrh9fjD6/kyj7yhV/GVNT1y9jetrvn5WTHW36rmE5+rFrMtPmazt89305dZyENn068bVtce2m797xxRjXPvr882nw4sXZabOSqfJPAbTkPV+Rs1WNBR5vfavzGhtn2nyrWH3bVZqaru+Xw7d31VdG3PSqz7/iRT+Ocn6S0fclGXF+APFTtz7/4q51isy5Nn3XJU4TPh37ftijoSHeHqd7X9k+/biFdaXgi9HndxJlL1zo/5/iWOPKEKm5H0ypmjZ917Ac/fjwtR3l3Pfltm3l6MenE1/bFtdu+u4dX4zlfAmVT4NLlpReRj6VJO73AXSISLuI1AN4AsCR6/xPyQwMVM9eK2VUs+xa8i8Fal6btVQ/89m/FJgX2qx1/VSz7LnQ5qxXlQOAiKwH8CsAdQB2q+rOa51f6spdki3meuVu7j2pTXJdqE1Sq1xLmxV9w7CqvgbgtUrKIKQaUJukVqE2SaWk++tghBBCCCkLJm5CCCEkIJi4CSGEkIBg4iaEEEICgombEEIICQgmbkIIISQgmLgJIYSQgGDiJoQQQgKCiZsQQggJCCZuQgghJCCYuAkhhJCAYOImhBBCAoKJmxBCCAkIJm5CCCEkIJi4CSGEkIBg4iaEEEICgombEEIICQgmbkIIISQgRFXn7s1E/gPgbMxLywH8d84cSYf5HOMPVLV5rp0pB2pz3sYYqjbnQ50B8yPOsrU5p4nbh4gcV9W70/ajmjDGMMliTMUwxvDIWjw+5kOcs4mRU+WEEEJIQDBxE0IIIQFRK4n7t2k7MAcwxjDJYkzFMMbwyFo8PuZDnGXHWBOfcRNCCCGkNGplxE0IIYSQEmDiJoQQQgIi1cQtIg+LyN9EZExEtqfpS5KIyG4RuSAip/Jsy0RkVEQ+zP1dmqaPlSAiLSLyJxE5LSJ/FZGf5+xZipHaDBBqM0yyrksgWW2mlrhFpA7AbwB0A1gJoFdEVqblT8IMAni4yLYdwBuq2gHgjdzzUJkEsE1VVwJYDeBnubrLRIzUZpj1loPaDJNBZFuXQILaTHPEvQrAmKp+pKpXALwMYGOK/iSGqr4NYLzIvBHAUO7xEIAfzalTCaKq51T1z7nHEwA+AHATshMjtRko1GaYZF2XQLLaTDNx3wTgn3nP/5WzZZUbVfVc7vG/AdyYpjNJISJtAH4I4BiyEyO1mQGozeDJSp3NoFJtcnFaCqjtwQt+H56INAEYBvALVb2U/1pWYpxvZKXeqM1skaU6S0KbaSbuTwG05D1fkbNllfMi8n0AyP29kLI/FSEiC2Hi26eqf8iZsxIjtRkw1GZmyEqd/Z+ktJlm4n4fQIeItItIPYAnABxJ0Z9qcwTA5tzjzQAOp+hLRYiIAPg9gA9UdVfeS1mJkdoMFGozU2SlzgAkrE1VTe0AsB7A3wH8A8Av0/Ql4bgOADgH4CrsM6h+AN+DrRj8EMAfASxL288K4lsLm875C4CTuWN9xmKkNgM8qM0wj6zrMhdjYtrkV54SQgghAcHFaYQQQkhAMHETQgghAcHETQghhAQEEzchhBASEEzchBBCSEAwcRNCCCEBwcRNCCGEBMT/AFJG3W6q8B5gAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 504x216 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "G1 = ot.emd(a, b, M1)\n",
    "G2 = ot.emd(a, b, M2)\n",
    "Gp = ot.emd(a, b, Mp)\n",
    "\n",
    "# OT matrices\n",
    "pl.figure(3, figsize=(7, 3))\n",
    "\n",
    "pl.subplot(1, 3, 1)\n",
    "ot.plot.plot2D_samples_mat(xs, xt, G1, c=[.5, .5, 1])\n",
    "pl.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')\n",
    "pl.plot(xt[:, 0], xt[:, 1], 'xr', label='Target samples')\n",
    "pl.axis('equal')\n",
    "# pl.legend(loc=0)\n",
    "pl.title('OT Euclidean')\n",
    "\n",
    "pl.subplot(1, 3, 2)\n",
    "ot.plot.plot2D_samples_mat(xs, xt, G2, c=[.5, .5, 1])\n",
    "pl.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')\n",
    "pl.plot(xt[:, 0], xt[:, 1], 'xr', label='Target samples')\n",
    "pl.axis('equal')\n",
    "# pl.legend(loc=0)\n",
    "pl.title('OT squared Euclidean')\n",
    "\n",
    "pl.subplot(1, 3, 3)\n",
    "ot.plot.plot2D_samples_mat(xs, xt, Gp, c=[.5, .5, 1])\n",
    "pl.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')\n",
    "pl.plot(xt[:, 0], xt[:, 1], 'xr', label='Target samples')\n",
    "pl.axis('equal')\n",
    "# pl.legend(loc=0)\n",
    "pl.title('OT sqrt Euclidean')\n",
    "pl.tight_layout()\n",
    "\n",
    "pl.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Dataset 2 : Partial circle\n",
    "--------------------------\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAbQAAADSCAYAAAAi7jW1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAYzUlEQVR4nO3df7QcZX3H8fen/FJAC+WmKJAQFUrNUVrlinhqT71A24AcsFQjtHpAoSAt1bbUHCwmuQm1CrRUrTYW8QeCNUbEir8K6A3+KmAuFRRENFAQIj9CLL+ECpFv/5jnhmXZ2d17d3ZndvbzOmdOdveZO/Ps7Ob57vPM95lRRGBmZjbsfqXsCpiZmRXBAc3MzGrBAc3MzGrBAc3MzGrBAc3MzGrBAc3MzGrBAc1sliQdL+lbZdejSJJC0j7p8YckLStouwskPSxpm/T8SkknFrHttL2vSDquqO3ZcHNAs1mT9EpJ/yXpAUk/k/RtSS8ru15VIGlhCg7blliH2yQdOte/j4i3RMSZRewnIn4SETtHxC/nWp+G/U1Kuqhp+4dFxAW9btvqobT/dDacJD0b+CJwCrAW2B74XeAXfdjXthGxpejtlq2u76vZqLxPqw730Gy2fgMgIj4VEb+MiEcj4vKI+B6ApF+R9E5Jt0u6V9InJP1qKnuVpDsbN9b4Kz/9Ar9Y0kWSHgSOl7SNpL+TdIukhyRdK2l+Wv83JV2Reok3S1qSV2lJb5J0U9rGrZJObih7laQ7JZ2W6nyXpDc1lO8m6VJJD0r6DvCCNsfnG+nf+9NQ2yvSEOW3Jf2zpM3ApKQXSJqStFnSfZI+KWmXhn2+VNJ3U30/I+nTkv6+ofwISddJuj/1lvdPr18ILAC+kPa/NOd4vD29z59KenNT2cdn9iVpTNIX035+Jumb6TN+2n4aeqcnSPoJMJXTY32BpO+k4/l5Sb/W+Dk01eU2SYdKWgz8HfD6tL/rU/nWIcwO372Zehwn6SfpmJ/RsJ8DJU2nOt0j6dw2n7FVVUR48dL1Ajwb2AxcABwG7NpU/mZgA/B8YGfgEuDCVPYq4M6m9W8DDk2PJ4HHgdeQ/dh6JvB24PvAfoCA3wJ2A3YC7gDeRDbS8BLgPmBRTr1fTRaIBPwe8Ajw0oZ6bQFWAdsBh6fyXVP5GrLe6E7Ai4CNwLdy9rMQCGDbhteOT9v/y1TXZwL7AL8P7ADMIwuE703rbw/cDrwt1edo4DHg71P5S4B7gZcD2wDHpeO4Q/MxzanjYuCe9F52Av491XmfVP7xhn29G/hQqsd2ZL1xtdpPw3v/RNruM5uPB3BlOn4z+/4scNEsvh8XNZVfCZzYxXdvph4fTvX6LbJRhRem8quAN6bHOwMHlf1/zcvsF/fQbFYi4kHglTzZOGxKvZfd0yp/CpwbEbdGxMPAO4Bj1P05pasi4j8i4omIeBQ4EXhnRNwcmesjYjNwBHBbRHwsIrZExHfJGsfX5dT7SxFxS9rG14HLyRrnGY8DqyLi8Yj4MvAwsJ+yZIY/BpZHxM8j4gayYD5bP42If0l1fTQiNkTEFRHxi4jYBJxLFmgBDiILfO9P9bkE+E7Dtk4C/i0iromsl3wBWeN8UJd1WQJ8LCJuiIifkwWKPI8DzwX2TnX5ZkR0ugDsZDpWj+aUX9iw72XAknSce9XNd29lOv7XA9eTBTbI3uc+ksYi4uGIuLqA+tiAOaDZrEXETRFxfETsRfZLew/gval4D7LexYzbyRrn3enOHU3P5wO3tFhvb+DlaSjsfkn3kzVoz2m1UUmHSbo6DZvdT9YLG2tYZXM89XzPI2S/1Oel+jfWq/H9desp70vS7pLWSNqYhlcvaqjPHsDGpsDR+Pd7A6c1vff56e+6sQfdv59zyHo9l6eh2tO72H7zZ9iu/Haynt9Yzrqz0c137+6GxzOfMcAJZMPpP5S0XtIRBdTHBswBzXoSET8kG6J6UXrpp2QN7owFZMNt9wA/B3acKUi/yuc1b7Lp+R20Pmd1B/D1iNilYdk5Ik5pXlHSDmS9t38Edo+IXYAvkw0/drIp1X9+03vKk9d7aX79H9JrL46IZwNvaKjPXcCekhrr17j/O4B3Nb33HSPiUx3qMOMuunw/EfFQRJwWEc8HjgT+RtIhHfbTaf/N+36cbLi40/ej03bbfffaiogfR8SxwK8DZwEXS9qp099ZtTig2awoS8Q4TdJe6fl84FhgZojmU8BfS3qepJ3JGu5Pp97Pj4BnSHq1pO2Ad5KdQ2rnfOBMSfsqs7+k3cgyLX9D0hslbZeWl0l6YYttbJ/2swnYIukw4A+6eb+RpZtfQpbIsaOkRWTnrPJsAp4gO4/TzrPIhjUfkLQn2bnCGVcBvwROlbStpKOAAxvKPwy8RdLL0zHZKR3TZ6Xyezrsfy1Zws0iSTsCK/JWVJZ8sk8Krg+kej3R5X7yvKFh36uAi9Nx7vT9uAdYKCmv3Wr33WtL0hskzYuIJ4D708tPtPsbqx4HNJuth8iSEa6R9HOyQHYDcFoq/yhwIVmSw/8A/0eWDEFEPAD8OVmQ2kj2i/wpWW0tnEvWAF8OPAh8BHhmRDxEFpSOIftlfjfZL+unBci07lvTdv4X+BPg0lm851PJhqbuJuuNfixvxYh4BHgX8O00HJh3Xmsl8FKyIPElsqA5s43HyBJBTiBrXN9AFsB/kcqngT8DPpDezwayxJMZ7wbemfb/ty3q+BWyIeKp9LdTbd77vsBXyYLvVcC/RsS6bvbTxoVkx/Fu4Blkn00334/PpH83S/rvFtvN/e51YTFwo6SHgfcBx7Q5B2gVNZOtZGYVJuka4EMRkRtMzUade2hmFSTp9yQ9Jw05HgfsD/xn2fUyqzJfKcSsmvbjyblvtwKvjYi7yq2SWbV5yNHMzGrBQ45mZlYLDmhmZlYLlT2HNjY2FgsXLiy7GmZmViHXXnvtfRHRfEEGoMIBbeHChUxPT5ddDTMzqxBJuZdq85CjmZnVggOamZnVQiEBTdJH0w31bsgpl6T3S9og6XuSXlrEfs1q4+yzYV12RanJyfTaunXZ67NZx2yEFdVD+zjZtdDyHEZ2Tbh9ye7ltLqg/ZoNh07B6GUvgyVLYN06Vq5MZUuWZK/P6LSOA56NuqLuFEp2R9gbcsr+DTi24fnNwHPbbe+AAw4Is9qYmooYG4uYmgp46vPmdVay7Oll3azTzT7MhhwwHXlxKK9gtkuHgPZF4JUNz78GjLdY7yRgGphesGBBnw+LWcHOOmtr8FixIr02NZW9PvM4JxitWJH9b1zJsgiIlSwLaNhOl+vk7qNT3cyGxNAEtMbFPTQbOm16SD0Foxb7mHVQdO/NaqIKAc1DjjYauhgSnPNwYa/Dlu3K3IOzIVGFgPZq4Ctkt5g/CPhOp+05oFnldGj0e+ohdRNQuhzSnFMP0T04GxJ9D2hktz6/C3ic7A6zJwBvAd6SygV8ELgF+H6n4cZwQLMq6qWHNIgeUA/n8LoqN6uAgfTQil4c0KyShjXLsEPdujrHZ1YB7QKarxRi1qjNXK7JSdDBE6y67xSWcyar7jsFHTzx5Hrr18PatTAxwYoVwMRE9nz9+oG/jafpULfJSYipdSwfW80qlrF8bDUxtS57b57fZsMiL9KVvbiHZqXo1Muq67Bcu/dd5Z6njRw85Gg2C3lBq84Ne6/n38wGxAHNrEttzyWNaGq7z69ZlTigmTVyb2T2fAUSqwgHNLNGPl80Oz5eViHtApqzHG30zGT4LVnCSpZnV6xPGYCVzlQsS7tj0u5Ymg1aXqQre3EPzfrF54SK42Npg0abHpqy8uoZHx+P6enpsqthdZXuJbbqvlNYPrbavYpe+FjaAEm6NiLGW5V5yNHqqd1k4JkbY65dywpWbR0ym1nfZqHdsfSEbBswBzSrp3Z3d/Z5suK0O5bd3IXbrEAecrT68lBY+fwZWME85Ggjp+N1F63v/BnYoLmHZvXl3kH5/BlYwdxDs9HjxI/y5X0GJ5/sZBHrCwc0G36tsunWrIGjj3biR5nyEkbAySLWFx5ytOHX0BPQwRPE1DpfsaLqPBRpc+QhR6s3X35pqDhZxPrFAc2GnhvI4dL27thmPXBAs6HnBnLIOGHH+sQBzYafG8jh0ipZ5Oijs0QenPloc+eAZsMj79qA55zjS1kNk6VLt57f3Po5HnMMXHKJMx+tJ85ytOHhbMZ6c+ajdcFZjlYPzmasLSf2WBEc0GxouNGrLyf2WBEc0GxouNGrMSf2WAEc0Gx4uNGrL9+jzgrggGbDw41efTVnPs6k6y9d6jR+61ohAU3SYkk3S9og6fQW5cdL2iTpurScWMR+reaa0/SXLt36+tZGbmLiydetPny3a5uDngOapG2ADwKHAYuAYyUtarHqpyPit9Nyfq/7tRHgRm10OaPV5qCIHtqBwIaIuDUiHgPWAEcVsF0bdW7URpYzWm0uighoewJ3NDy/M73W7I8lfU/SxZLmt9qQpJMkTUua3rRpUwFVs2HmRm10OaPV5mJQSSFfABZGxP7AFcAFrVaKiPMiYjwixufNmzegqllVuVEbYc5otTkoIqBtBBp7XHul17aKiM0R8Yv09HzggAL2a3XnRm10OaPV5qCIgLYe2FfS8yRtDxwDXNq4gqTnNjw9EripgP1a3blRG11O47c56DmgRcQW4FTgMrJAtTYibpS0StKRabW3SrpR0vXAW4Hje92v1ZDT9C2PM16tC77avlWHr6Zv7fhq/Iavtm/Dwmn6lsMZr9YNBzSrDDdalscZr9YNBzSrDDdalssZr9YFBzSrDjdalscZr9YFBzSrDjdalqcxjX/HlA07McHkIynj1Sn8hgOala0xVf+R1GitW5c1WuA0fXs6p/BbDgc0K5cbJ5stZ8NaDgc0K5cbJ5slZ8NaHgc0K5UbJ5stZ8NaHgc0K5UbJ5s1Z8NaDgc0K5cbJ5stZ8NaDl/L0cp19tlZAshENsw4OUkWzNavd3ajmT1Nu2s5OqCZ2XDyj6GR5IsTW/U03yoGPDnWZsdTPqyJA5qVw42R9cpTPqyJA5qVw42R9chTPqyZA5qVwo2R9cpTPqyZA5qVwo2R9cxTPqyJA5qVw42R9crz0ayJ0/atHE65NrM58Dw0MxsJW38cWW15HppVi+egWZ+sXFl2DaxMDmg2eJ6DZmZ94IBmg+c5aFagyUmQsgWefOyhx9HjgGYD5zloVqTJSYjIFnjysb9Po8cBzQbOc9DMrB8c0GzwPAfN+mTFirJrYGVyQLPB84RY6xP38kdbIfPQJC0G3gdsA5wfEe9pKt8B+ARwALAZeH1E3NZum56HZmZmzfo6D03SNsAHgcOARcCxkhY1rXYC8L8RsQ/wz8BZve7XhpjnoZlZHxQx5HggsCEibo2Ix4A1wFFN6xwFXJAeXwwcIs0k2drI8Tw0M+uDIgLansAdDc/vTK+1XCcitgAPALs1b0jSSZKmJU1v2rSpgKpZJXkempn1QaWSQiLivIgYj4jxefPmlV0d6xPPQzOzfigioG0E5jc83yu91nIdSdsCv0qWHGIjyPPQzKwfigho64F9JT1P0vbAMcClTetcChyXHr8WmIqqXubf+s/z0MysD3oOaOmc2KnAZcBNwNqIuFHSKklHptU+AuwmaQPwN8Dpve7XhpjnoZlZH/h+aGZmNjR8PzSrFs9DM7M+cECzwfM8NDPrg23LroCNoKfMQzsFlqz2PDQz65l7aDZwnodmZv3ggGYD53loZtYPDmg2eJ6HZmZ94IBmg+d5aGbWB56HZmZmQ8Pz0KyaPB/NzArkgGbl8Xw0MyuQ56FZeTwfzcwK5B6alcbz0cysSA5oVhrPRzOzIjmgWXk8H83MCuSAZuXxfDQzK5DnoZmZ2dDwPDSrNs9HM7MCOKBZ+TwfzcwK4HloVj7PRzOzAriHZqXzfDQzK4IDmpXO89HMrAgOaFY+z0czswI4oFn5PB/NzArggGblW7p0awLI5I4phX9igslHlmblTuE3sy44oFm1OIXfzObIaftWLU7hN7M5cg/NKsUp/GY2Vw5oVilO4TezueopoEn6NUlXSPpx+nfXnPV+Kem6tFzayz6t5pzCb2Zz1GsP7XTgaxGxL/C19LyVRyPit9NyZI/7tDpzCr+ZzVGvAe0o4IL0+ALgNT1uz0ZdYwr/JE+m6y9d6ivxm1lbvQa03SPirvT4bmD3nPWeIWla0tWScoOepJPSetObNm3qsWpWC07jN7MudUzbl/RV4Dktis5ofBIRISnvbqF7R8RGSc8HpiR9PyJuaV4pIs4DzoPsBp8da2/15zR+M+tSxx5aRBwaES9qsXweuEfScwHSv/fmbGNj+vdW4ErgJYW9A6s1p/GbWbd6HXK8FDguPT4O+HzzCpJ2lbRDejwG/A7wgx73ayPCafxm1q1eA9p7gN+X9GPg0PQcSeOSzk/rvBCYlnQ9sA54T0Q4oFl3nMZvZl3q6dJXEbEZOKTF69PAienxfwEv7mU/NsLapfH7PJqZNfCVQqzamtP4IQtmKcvRqfxmNsMBzYaPU/nNrAVfbd+Gj1P5zawF99Bs6DiV38xacUCzoeNU/ho4++ynZ6r6PKj1yAHNho9T+Ydfw3lQwOdBrRAOaDZ8WqXyH300rFkDOPNxKDScB2X58q0/UHwe1HqhiGpeMnF8fDymp6fLroYNi4Zemw6eIKbWuZEcBsuXw5lnwrJlsGpV2bWxISDp2ogYb1XmHprVw1MyH/2LfyisWwerV2fBbPVqDxlbzxzQrBac+ThkGnrUrPJ5UCuGA5rVgjMfK6ohm/Ep5zbPOeepPWjfmdwK4HNoVg8+h1ZN/lysYD6HZvWXdxHjc85p3UNw9uNg+NymDZADmtVDq4sYT0zA29/u6z6WyOc2bZA85Gj1l4LYqvtOYfmYr/s4cD7+ViAPOdrIcg9hQPKSP04+2Vd1sYFxQLNay81+3DGnAfa5tbnJu6UP5N+g1axoEVHJ5YADDgiznk1NRYyNRUxNBTQ8/6d/av361FTZNR5e6RiuZJmPpfUNMB05ccM9NKu3vOzHLVucfVcgD+1aJeRFurIX99Csn1asiIDIehMQK1kWkL1uOc46a2uva+txmprKXp957B6a9RltemilB668xQHN+i6vAe7UcI+qvOHbqan2ZWYFahfQPORoo6ndPdXyEhxGfe5au0nSeUO7Tv6wQcqLdGUv7qFZX3n4bNY8TGtVgIcczbo30g23A71VnAOa2Wy1a7jrfI7N58ms4toFNJ9DM2vW7vwaDPc5trwresxMKPd5MhtmeZGu7MU9NCtNNz2wYc2Q7NDLGunhVhsKeMjRrDhtG/2yh+V6CcbdlpuVqG8BDXgdcCPwBDDeZr3FwM3ABuD0brbtgGaV1q7R7+X8W6/lvfbAyg7IZh30M6C9ENgPuDIvoAHbALcAzwe2B64HFnXatgOaVVabRr/ngNFrecM6c+qBVX3I1EZe34ccOwS0VwCXNTx/B/COTtt0QLPK6jW1vY/l7oFZ3ZUd0F4LnN/w/I3AB3LWPQmYBqYXLFjQ58Ni1gc9Dvn1Wt5YB/fArI56CmjAV4EbWixHNaxTSEBrXNxDs6E0iKSMLs7fuQdmdVV2D81DjmYz+n0OzT0wq7l2AW0QE6vXA/tKep6k7YFjgEsHsF+z6uk0ObnX8qVLt97TbevE6YmJ7HWzmlMW8Ob4x9IfAf8CzAPuB66LiD+UtAfZMOPhab3DgfeSZTx+NCLe1Wnb4+PjMT09Pee6mZlZ/Ui6NiLGW5Vt28uGI+JzwOdavP5T4PCG518GvtzLvszMzNrxtRzNzKwWHNDMzKwWejqH1k+SNgG3l12PHowB95VdiSHg49QdH6fu+Dh1Z5iP094RMa9VQWUD2rCTNJ134tKe5OPUHR+n7vg4daeux8lDjmZmVgsOaGZmVgsOaP1zXtkVGBI+Tt3xceqOj1N3anmcfA7NzMxqwT00MzOrBQe0PpF0jqQfSvqepM9J2qXsOlWVpNdJulHSE5Jql3nVC0mLJd0saYOk08uuT1VJ+qikeyXdUHZdqkzSfEnrJP0g/Z97W9l1KpIDWv9cAbwoIvYHfkR2lwFr7QbgaOAbZVekSiRtA3wQOAxYBBwraVG5taqsjwOLy67EENgCnBYRi4CDgL+o03fKAa1PIuLyiNiSnl4N7FVmfaosIm6KiJvLrkcFHQhsiIhbI+IxYA1wVMl1qqSI+Abws7LrUXURcVdE/Hd6/BBwE7BnubUqjgPaYLwZ+ErZlbChsydwR8PzO6lR42PlkrQQeAlwTbk1KU5PV9sfdZK+CjynRdEZEfH5tM4ZZN38Tw6yblXTzbEys8GQtDPwWeCvIuLBsutTFAe0HkTEoe3KJR0PHAEcEiM+P6LTsbKWNgLzG57vlV4zmzNJ25EFs09GxCVl16dIHnLsE0mLgaXAkRHxSNn1saHku71boSQJ+AhwU0ScW3Z9iuaA1j8fAJ4FXCHpOkkfKrtCVSXpjyTdCbwC+JKky8quUxWkpKJTgcvITt6vjYgby61VNUn6FHAVsJ+kOyWdUHadKup3gDcCB6d26TpJh3f6o2HhK4WYmVktuIdmZma14IBmZma14IBmZma14IBmZma14IBmZma14IBmZma14IBmZma14IBmZma18P9MORVkuis1FwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 504x216 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfAAAACsCAYAAACEl/7eAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO2de5RkVX3vP7+q7ume6emZnhnICDMoKvgA4vsBPgnqFRCFlRgjGgKGBM0yV1ni+y5vTDQGcxPRXE3MKAZ8Ar4R9Ublcb0mioCiRlBBBHkMj2GmoZlnd9Xv/rH36dpVU11V3fU81d/PWr361D6vffb5nr3P73e+dcrcHSGEEELki0K/KyCEEEKIxaMBXAghhMghGsCFEEKIHKIBXAghhMghGsCFEEKIHKIBXAghhMghGsBrMLNjzeyO5PPPzezYVpYVopuYmZvZYT3c36FxnyPx8zfN7PRWlhViKZjZu83s03H64Wb2kJkVmy27XMn1AG5mt5rZ7niSs78Pd3If7n6ku1/VyW0OGmZ2hpl9r9/16AVm9hwz+08ze8DMtpvZf5jZ0/tdr3Yxs6vMbE/NtfC1Tu7D3U9w9ws7uc1BI6835Z3WddTTnzWYn92wPVTz90dL3Wct7v5bd1/t7qVObXMQMbMLzOy9S1l3GO6WX+ru3+l3JcTgY2ZrgMuAvwAuAVYAzwX29qEuxS50TH/p7h/v8DbFgNNJXZuZAbaIVabcfW6x+xGdIdcReCNq0yt10oHrzezfzOwuM9thZl9ZYDu3mtkL4/TKeLe0w8xuAJ5es+zBZvZFM7vPzH5jZm9I5j3DzL5vZtNmttXMPmxmK5L5bmavM7Ob4jIfiRdTvToVzeydZvZrM5sxs+vM7JA471lmdk28E7/GzJ6VrHeGmd0S1/mNmb3azB4PfBQ4Jt5BTy+hufPCYwDc/XPuXnL33e7+LXf/Kcy36z+Y2bbYTq+v0cy8FuLnWo193szujm3/XTM7Mpl3gZn9i5l9w8x2Ar9nZmNxf781s3vM7KNmtjJZ5y1RK3eZ2Z8u9aDrZVgsScdHXf+jmd0W6/69tB7JOvNRWW1bAS+pWXatmZ0f63+nmb3XYirUzB5tZleY2f1x/c+Y2VSy7q1m9mYz+2msz8VmNt7g+P7czG6Mur7BzJ4Syx8f6zxt4VHYy5J1TozLzsT6vdnMJoBvAgdbJaI8eAlN3mva1fVVZva3ZvYfwC7gU4QbgA/bErOaVhPB12rQzI40s29byBbcY2bvrLON2j77kWb2f+M5+zZwQM3yR1vIQkyb2U8sefRpZq9JNHKLmb02mXesmd1hZueY2b1Rs69pcGwLjh1RizfH47o0048Fzovbf9DMfmZmR5nZWcCrgbfaUrJm7p7bP+BW4IULzHs38Onk86GAAyPx89eBi4F1wCjw/Fh+LHBHvX0A5wL/D1gPHAL8V7Ys4WboOuB/Eu6AHwXcArw4zn8qcDQh63EocCNwdrIfJ9xFTwEPB+4Djl/g2N4C/Ax4LOFu+YnAhlivHcBpcT+nxs8bgAngQeCxcRsHAUfG6TOA7/X7fPZAL2uA+4ELgROAdTXzXwf8Ip7b9cCVNZqp0lsdjf0pMAmMAR8Erk/mXQA8ADw7amUcOA+4NO5rEvga8Hdx+eOBe4Cj4rn7bKzLYQsc21XAny0wb7/zm24L+EhcfxNQBJ4Vj+HQmuOf30cLbfVl4F9j3X8H+CHw2jjvMOBFcR8HAt8FPlhzzf0QODhu+0bgdQsc2x8CdxJupi1u+xGEa/pm4J2E6/E4YCbR/1bguXF6HfCUetd/Hv5oX9dXAb8FjiT0G6ON9BTXqdJGMz2mGiRofStwDuE6mASeWXtN1dHf94EPRN08L57PbNlNsQ1OJFxfL4qfD4zzXwI8Omrk+YQblfSczwF/E4/9xDh/3QLHttDYcRywDXhKrOP/Br4b572YMD5MxTo8Hjgo6Rveu6Rz32/xtSncW4GHgOnk789rhVArBsLgVa53gmg8gN9CMqgCZ1EZwJ8J/LZmW+8A/m2Bup8NfDn57MBzks+XAG9fYN1fAifXKT8N+GFN2fcJF89EbJ8/AFbWLHMGy2AAj8f6+HjB3BEv2kuBjXHeFSQDBfDfWMQAXrOfqbju2vj5AuCTyXwDdgKPTsqOAX4Tpz8BnJvMewzNB/BdNdfCexY6v9m2CJ3dbuCJdbY5f80k+8gG8AXbCthISN+uTOafCly5QN1PAX6cfL4V+OPk898DH11g3X8H3lin/LnA3UAhKfsc8O44/VvgtcCamvWOJWcDeAd0fRXwN3X01MoAPl3z9/h661M9gJ+anu+a7b6bOgM4IaiZAyaSZT+bLPs24FN1tHH6Avv5SqabeM53k9yMAPcCR9dZr9HYcT7w98nn1cBsPI7jgF8RgrhCzXoXsMQBfBhS6Ke4+1Ty97EW1jkE2O7uOxa5r4OB25PPtyXTjyCk3qazP8Ld/0YAM3uMmV1mIcX6IPA+alJAhA4nYxdBAAvV/9cL1O+2mrLbgE3uvhP4I8Ld+FYz+7qZPW6hAx1W3P1Gdz/D3TcTotuDCdEyND6/DYlpynMtPNZ4kDAIQfU5Trd9ILAKuC7Ry/+J5UutyxtqroV3tbDOAYQoqJ6eGtHsWhgl6Cw7tn8lROKY2UYzuyimrh8EPk13roXb3b1cU8dNcfoPCJHWbTEte8wC288FHdD17XXKWuGAGs3d2MI6C52zRhwM7Ij9WEat5v6wpv99DmHAxcxOMLMfxNT2NOHcp5q736uf5S+kuUZjR1X/6+4PEbIAm9z9CuDDhGzXvWa2xYJ3oS2GYQBfiJ2EDjLjYcn07cB6S567tchWwgnMeHjNNn9TI+ZJdz8xzv8XQhrrcHdfQxjcF2MWSbmdkA6q5S6CkFMeTkgx4u7/7u4vIoj6F0B2s+NLrEeucfdfEO5+j4pFjc4vNNbUq4CTgRcCawl33VB9jtN23ka46z8y0ctad886jWZ1WQxV9TaztN7bgD3U11Mjml0Le6nu3Ne4e+YJeB+hLX43Xgt/THeuhUPMLO3j0mvhGnc/mXBT8RVCxguG4FpYgq5h/+Nutx2a9b+PWuT2tgLrok8ho1Zzn6rpfyfc/VwzGwO+CPwDISsxBXyDpWmu0dhR1f/Gum6gorl/cvenAkcQMmpviYsuua2HeQC/Hniehe8SriWkswFw960Es8o/m9k6Mxs1s+e1sM1LgHfEdTYD/z2Z90NgxszeZsEUVIwmhczoNkl4Bv1QjHz/oo1j+zjwHjM7PJojnmBmGwiifIyZvcrMRix8peMI4LIY9ZwcRbWX8Oghi07uATZbYqobRszscdGosjl+PoSQzvtBXOQS4A1mttnM1gFvr9nE9cAro16eBrw8mTdJaNf7CR3X+xrVJUaGHwPOM7MsMt1kZi9O6nKGmR1hZquAv1raUQPwE+BIM3uSBTPYu2vq8QngAxZMmEUzOyZ2eo1YsK3i9fUt4B/NbI2ZFSwY154fF5kk6O8BM9tEpSNbCh8H3mxmT43XwmFm9gjgakIU9dZ4vo4FXgpcZGYrLBg417r7LOG6TK+FDbHPyAUd0HU97mHxg2zK9cDvm9kqC2bJM5N5lwEHmdnZFoyck2b2zEYbc/fbgGuBv47n7zmE85nxaeClZvbiqOFxC+a0zQQPxBjBVzRnZicQHiMsmiZjx+eA18TrbIzQB1zt7rea2dPN7JlmNkq4udlDteaW1NbDMIB/zaq/h/hlAHf/NsFo8FOCeeCymvVOIzyf+AXhecfZLezrrwkpkt8QOqhPZTM8fCXoJOBJcf42QueSdQRvJkRpM4SO++JFH2mFDxAuym8ROp/zCc8b7491OIcwkLwVOMndtxHO9ZsId4nbCUaO7CbiCuDnwN1mtq2Neg06MwSvwtUWnOA/IBgRz4nzP0Z4bvYT4EfAl2rWfxch2ttB0MJnk3mfJGjjTuAGKp1nI95GMFr9IKaSv0MwJuLu3ySkQK+Iy1zRwvY+XHMtXBe39SuCQec7wE1A7Xf+30wwRV5D0Mb7ad43NGurPyF0nDcQ2usLxHQmoe2eQjD1fb3Oui3j7p8H/pZwLmYI0fR6d99H6OBPIFyL/wz8SYxOIVz/t8Z2fx3BCZxFr58Dbomp2Dy40NvVdT0+BLzcgsv6nxosN12juTfF8vOAfYTB6ULgM9kK7j5DMJm9lPCo5Cbg91qo06vicW4n3NB+Mtnm7YQM2DsJA/XthBvDQtzfGwh95o64nUtb2N9C1B07PHyd+V2EaH8roa94ZVxnDeE87CD0E/cD/yvOOx84Iuqt7rehFsLiQ3QhRA1mdijhZmzU9V1XMSRI18PDMETgQgghxLJDA7gQQgiRQ9oawM3seDP7pYU3z7RijBCip7SjUXe/1d1NaUbRLfrRh0rXw8OSn4FbeC3irwhmhDsIBphT3f2GzlVPiKUjjYpBRvoU7dJOBP4M4GZ3vyU6Pi8iuACFGBSkUTHISJ+iLdr5NbJNVL+95w6CxX9BVhTGfWVxEoqVn3f1kTDtxcp36rPpcvIrsF7I5iUbLFTPA6AYMgpWqGQWCoXwdbtiUjZioWykUNq/zJIy9i8rxuWKyffvC3G6mPz+iM3/b1xWD4/b8zpl6aulyjGDUkq2V4oNU0oaZi5Oz1FMysL0XLnOcklZKU6Xy0mds+mkLDYLSVOx5+47trn7gfSHRWl0hY35OBNVr3eYP0+FRGSFemVhukrHcblqbYf/SfPOl9XT8fx/oBD1O1KsKGA06nc0afQV82WVDOkooSzTOFT0m2q2EI+3mT7L81qs1C/T4lyybqax2eTCzab3eaX7mY0X+2w93ZWShinF+iUXQXbodcsSLVo59g2lSp1ndt6VG30CjK6Y8PFV6yiNJpqKzeh1+stUP1mfmPWHAEWLGkjKMo0UU63UKyMr8/3KCslyWd+YyjsrsyrtxXrW6S9TmmmzFk80Wi/fXK7X13qmb0uWY7+yrI8tJUdXrlOW9atpn5xNz6VlsT8tJdfBnl9vXVCjXf85UQu/tnIWwHhhNcdM/T6sr7zEprQuvFhn37rKeyP2rQlK3LumchCzk+HAZpP38MxNxg5joiIWn4id16p982WrV4Vf1Vu7cs982Ybx8Ea+DWOVN/MduOIhAA4YnZkvW18MZRtGHpovmyrsiv8rv9a3KvYUk4XKyR2Pv0M/mvwe/UgcQIvWOPlRim+AnKPSA+2Jj6z2Jm+H3Bk7pZmkM5wuhx9umi5VXoS0vRRe8HXf3OR82bbZML1tb+WNgffvDQ28fU9l3endYXu7dlXOUWnnKACFncmx7QzHPvJQpQ1+8b43tfw60n5QpU9W8Ux7ATZSacts2sYrx24rww91+UTlB7t8IrRRaXVludnJ0Eb7JivnOpveN1lpo9nV2fKVLmRuTTjvhcnZ+bKJyaDfDRO75ss2rgpa3TRe+RG5TWNh+uDRytseHzbyAAAHFit6X1vINFs5h6viu3xSzdZjr4d67SpX6jcTdbm9NDpfdm/U3d1zlfeibJ1dB8Cdeyv9wF27w/y7d1beLrljZ2jf3TPJD5HNhPMxMlNp09GZ0JYrKpctozOhLcdmkpudh+JNzEylzpf/57sGWp9QrdGxlVM8+blv4KGDKxrdsyEc/741iX4mw3Fn/SHAionQJ06srPRba8bD9NTY7vmydWNBX2tHk7KRWDZS0d5UMUxPFir96mQhrLMmKVtVCO09ntxNrYqD/opkAB+LfeJoEmRkOiwkg2azvrOWUtJfZv1pKXl0vDf2q7PJEL4nzt/jlf3uKoc23+kVfc+Ug0YfLFc0OlMKZQ+k/e9c6FenZ5N+dTYsN72v0o9M7wnTM3sq/ch/nfyeBTXazgB+J9Wv59scy6pw9y3AFoA1tt5L928n7Rqy6epXgNV7CVTl/qyC1cwLb7sHmE22WBl6O0S9VssG83Jyy59E9/NkVU5u9+oJcr4sWW58/tAT70l251xOyyoXT7fILuP0COdo3On3gaYardUngM/V8fbUadJ6cUD9FhitU5aec6v5X9lSWpP0JdCdIdtiqtnQya9Krp96g/mYxWNKDyMbzIuzSWFnr77dhI5yruoibK1vGECW1IeOff0aeEn6S8ZZW6THn2XcKuyjzy9arDoVmebSuDdL4SVFcXaVBj2L/Fs7t1XLZbtLbibGsvZLPX1WuwJQmKuqZihrqQpdo53dXwMcbuE3WlcQ3jjTzttthOg00qgYZKRP0RZLjsDdfc7M/pLwir4i8Al3/3nHaiZEm0ijYpCRPkW7tPUM3N2/QfgBjZawkSLFqfWU7t8+X1as+Q9pOr1RKh3qpx73Tx1l6fSeptKhkk5vlEqH+SxNw1R6stx4Ve42HmliQplPp/cwlQ6VpNggpdIXq9H91u9pOr1eCjgxHcb/3Uulw/xZLFQ8JFk6vWEqHSqHkTwXr6TTu5NKhzSd3lrfMEgsVp++dhV7n/t0xr5+TaVwPp2edkgLP17seyodkmqlfWOWrk7diDWzSHToqcluken0NHMf0+ljaftl6fSqCzuuVEj6hKwKfZLW4ClaCCGEEE3pugu9imIR1k9VRSZZNC5jG0NpbMsVFlzn9aJuGdtkbBsESqMWHOiJiW0+GpexbdkZ2/KrZCGEEGIZowFcCCGEyCE9TaH7SJHSuom66XIZ2xhKY1ueMKzqJS5QP3UuY9swGdvyhY9kL25JdBpT5zK2seyMbYrAhRBCiBzS2wi8aOxbN1Z1r1cv2paxjaExtuWKQiG8MrVOs/Q/Ek+3JGNbM1o3tuULL2avTE2PIR6jjG3LztimCFwIIYTIIRrAhRBCiBzS+xT6miJp6jtL1sjYxlAb23JBweZ/aWyegUyny9jWKs2NbfnCC9kvjdXr32RsW27GtvwqWQghhFjG9DQCLxerf+M7EKJnGduQsa3fFAr4xMr6FqeBjMTTLcnY1oz6xracUXR8okS1yur1bzK2LQdjmyJwIYQQIodoABdCCCFySG9NbAWYnTTq3zfI2CZjW58pFPCJ8aqi/KTTZWxrldTYljes4KyY2FeV0q6oTMa25WZsUwQuhBBC5JCmEbiZfQI4CbjX3Y+KZeuBi4FDgVuBV7j7jmbb8iLMTkC9qLcaGdtkbGudTmnUi0Zp9Vjd6Dc/kXi6JRnbBoVOabRQKDOxcm9VWRYVy9gGy83Y1sriFwDH15S9Hbjc3Q8HLo+fhegXFyCNisHmAqRR0WGaDuDu/l1ge03xycCFcfpC4JQO10uIlpFGxaAjjYpusFQT20Z33xqn7wY2trRWAeYma1/Ev3/auoKMbTK2LZlFa9QLxuxkdXo63+l0GdsGnEVrtGjOmvG9defJ2AbLzdjWtonN3Z3q6ldhZmeZ2bVmdm1pZ+e7CiGa0UijqT5n90mfoj+0rNEHdvW4ZmKQWWoEfo+ZHeTuW83sIODehRZ09y3AFoCxhx/icxMLvcdXxrYqZGxrl5Y0mupz9fpDfN9kgXpRbb4j8XRLMrYNEIvW6NTjfsenxnY33bCMbTA0xrYGLDUCvxQ4PU6fDnx1idsRoltIo2LQkUZFWzQdwM3sc8D3gcea2R1mdiZwLvAiM7sJeGH8LERfkEbFoCONim7QNIXu7qcuMOsFi95byy/il7FtHhnbmtIpjXoRQgo9ZdjS6TK29YNOaXTEyqwba/05uIxtkHtjWwP0JjYhhBAih/T0XehWcEZX7asyfzW+K5SxrQoZ27pKuQD7FnxX/7BF4umWZGzLC0Urs3a0uYmtHjK2QS6NbQ1QBC6EEELkEA3gQgghRA7paQq9UCizetXequTWbMtpHRnb5pGxrSt4EWZXQ/NHN8OWTpexLS8Urcy6kfZe5iJjG+TL2LYwisCFEEKIHNLTCLxYcNaurA4jsvtjGduWyKAb23KEF2B2Ue/qH7ZIPN2SjG2DSNHKrG0zAk+RsQ0G3tjWAEXgQgghRA7RAC6EEELkkJ6m0EeszIbx+sk3GdvaZFCNbXmi6MytKVGtgFYf3QxbOl3GtkGkSJmpYud/kUzGNhhYY1sDFIELIYQQOaS3EXihxIax5vfsMra1ySAZ2/JE0SlMztbUPlOAjG3DaWzLF0VzJrv8NU0Z22CQjG2NUAQuhBBC5BAN4EIIIUQO6bmJ7cAVrSeMZWxrk4EwtuWHQsGZmNxTlQqunCcZ24bS2JYzipSZLCztx0wWi4xtMBDGtlaqKYQQQoj80DQCN7NDgE8CGwn3BVvc/UNmth64GDgUuBV4hbvvaLgzK3HA6MySKipjW5v0ydjWCzql0ZFimQ0T1V/RyWI/GdtgKI1tPaCTfWjByqzpw28NyNgGfTO2tVq9BZgDznH3I4Cjgdeb2RHA24HL3f1w4PL4WYh+II2KQUb6FF2h6QDu7lvd/Udxega4EdgEnAxcGBe7EDilW5UUohHSqBhkpE/RLRZlYjOzQ4EnA1cDG919a5x1NyE91GRnZdYX20sOy9jWJr02tvWYdjQ6WiixcVX9RzwyttUyJMa2HtNuH1rAWVXonwlPxjbotbGtpSo1w8xWA18Eznb3B9N57u4s8OTTzM4ys2vN7NqZHfl1f4rBZykaTfW5d0c+f8dc5INO9KHT23P6mmLRFVqKwM1slCC8z7j7l2LxPWZ2kLtvNbODgHvrrevuW4AtAIf/7krfMNK5mFLGtjbpgbGtVyxVo6k+Nx6x3jeNTzfdl4xtKfk2tvWKTvWhRz1hhY9baSC+PyRjG/TS2Na0KvUwMwPOB2509w8ksy4FTo/TpwNfXXIthGgDaVQMMtKn6BatRODPBk4DfmZm18eydwLnApeY2ZnAbcArulNFIZoijYpBRvoUXaHpAO7u32OBTB7wgsXsrGhlpgq7Ov7+Nxnb2qQHxrZu0imNriiU2DTWPIWeIWNbLfk1tnWTTvahBZxV5lQ/Llhy1TqCjG3QL2PbADxJEUIIIcRi6em70Is4U2lk14UayNjWJh02tuWJUZvj4NGGL8JaEBnbUvJjbMsbZsYKM6ovsizb0I8aVSNjG/TS2DYAp1wIIYQQi0UDuBBCCJFDeppCDwaMUrVBqks1kbGtTTplbMsRo5R42MgDbW1DxrZaBtvYljcKwJgVqDJIzadqZWzreyodOm5sa2lXQgghhMgPvTWxmTFZsBozVO+icRnblkg7xrYcMWJlDix2LgaUsS1lQI1tOcMwRinWnIwsWpOxbSiNba3uQgghhBD5QAO4EEIIkUN6mkI3YNyK1anWeTOUjG1DbWzLAUWctYUSnU7eythWyyAZ2/KFEdOrXlMIyNhWvb1hMba1tFkhhBBC5IceR+C2/8P5LFqTsW04jW05Ipgsi1TfJXcnGpexDQbD2JY/CrX9aBatydjGUBrbWt2cEEIIIfKBBnAhhBAih/TcxDay33cYIzK2DaexLUcUMFbZipo3dWW6lLGt2br5NrblA8PCD10kb+qqvL2rasGIjG25N7a1sgkhhBBC5Iem8ayZjQPfJYSCI8AX3P2vzOyRwEXABuA64DR3b/iS4crdY1Xh/sjYNjzGth7QKY1mJsuqn5ucj8ZlbJOxbWl0sg/NqPq5yRiNy9jGkBrbWlx1AfYCx7n7E4EnAceb2dHA+4Hz3P0wYAdwZqv1FaLDSKNikJE+RVdoOoB7IAvWRuOfA8cBX4jlFwKndKWGQjRBGhWDjPQpukVLSU8zKxJSPIcBHwF+DUy7e5aJuAPY1OpOq9M/2U7qLChj2/AY27pMJzWapiLn0+kytjUtk7FtYTrdh6bM96cytg2psa2FxRvh7iV3fxKwGXgG8LhW62NmZ5nZtWZ27X335/MFH2LwWapGpU/RC9SHim6wqFjJ3afN7ErgGGDKzEbiHeRm4M4F1tkCbAF42hPH93syX7l7TAplbBseY1uPWaxGm+kzi2RkbBuESDzdUj6Nbd3oQzNkbGNIjW1tLGZmB5rZVJxeCbwIuBG4Enh5XOx04KtLqacQ7SKNikFG+hTdopV49SDgwvgMpwBc4u6XmdkNwEVm9l7gx8D5XaynEI2QRsUgI32KrtB0AHf3nwJPrlN+C+FZTss4TsnL1ameiIxtQ25s6yKd0mgZZ6/PMmb7p3FlbKum/+n0/BjbOtmHtoqMbQyPsa0BA/CUQwghhBCLpadf+HFgjlLVHWDDaFzGNhnbekgZZ1d5tuoUNYrGZWwbhEg83VI+jW2LoVEWsx4ytjE8xrY6DMDpEEIIIcRi0QAuhBBC5JAep9CdPT7HuFUVAk1S6clyMrbJ2NYtyu7MeBnKs5XC2FwyttUwkOn0/BjblkqrjyHrIWMbQ2dsUwQuhBBC5JCeRuBlYK+XSe+X5qNxGdtqkLGt18xhbC+NQjGJwLNoXMa2HEXi6ZYaG9vySMkdLNFUgyxmPWRsY2iMbQPQ9EIIIYRYLBrAhRBCiBzS2xS6OzvLDoXUHBESEjK2pfNSZGzrFXNe5N7SaqqOKkuny9hWRX7S6Y2NbXkjvC1wjrG048rS6TK2LTtjmyJwIYQQIof0NAIvYcz4CJSTO/D5aFzGNhnb+susF7l7bm1NaTwqGdtyHomnW6pb61zgwCwOXmm7+WhcxrYqloOxbQCaWQghhBCLRQO4EEIIkUN6nEIvMF0eh0KSJ8vS6TK2IWNbf5n1Iltn1y0wV8a2lHyn03OcQndnj3v1IcR0uoxt9RlmY5sicCGEECKHtBxrmlkRuBa4091PMrNHAhcBG4DrgNPcfV+jbZS8wHRpVXVhFo3L2MYwGtt6RSf0uc9HuHPvVAt7k7Gtdgsp+YnEe0snNFrG2ONGVQc33zfK2NaMYTO2LaZJ3wjcmHx+P3Ceux8G7ADO7GTFhFgk0qcYdKRR0VFaGsDNbDPwEuDj8bMBxwFfiItcCJzSjQoK0QzpUww60qjoBq0miD8IvBWYjJ83ANPu8zmbO4BNzTYy5wW2l1bXnylj21Aa23pER/Q5Wy5y1+7a74E3Qsa2lHym03tGRzRaBnaVR6CQtk3s4GRsa5lhMbY1bUYzOwm4192vW8oOzOwsM7vWzK59aMds8xWEWASd1Ofe6d0drp0QndXoju3l5iuIZUMrceWzgZeZ2YnAOLAG+BAwZWYj8Q5yM3BnvZXdfQuwBWDzUWv9vrnJejgsqf4AAAWxSURBVItVI2Pb0BjbekDH9Ln6MQ/zu3euWWI1ZGyr3ULKMo/EO6bRxz1hzHf6aHWAOx+Ny9i2FPJsbGvafO7+Dnff7O6HAq8ErnD3VwNXAi+Pi50OfLXt2gixSKRPMehIo6JbtHP/8zbgTWZ2M+F5zvmdqZIQHUH6FIOONCraYlHJYHe/CrgqTt8CPGMx6895kW2zLaTQM2RsGwJjW+9oV5+lcoEdO1e2WQsZ21IGP53eW9rWqBeYKa+sbv6sS5SxrS3yaGwbgCcQQgghhFgsPX0X+ly5wLa9C3yNrBkytuXU2JYfyqUCu2fGO7hFGdtqt5CynCPxpVKiwIPlGo1mp0DGto6RF2PbADSVEEIIIRaLBnAhhBAih/Q2he4F7t870d5GZGzLmbEtR5QMZkbYTSfT6CBjWzWDlE7PG2UvMFNawGgpY1vHGXRjW057WiGEEGJ503MT2/Y9q5ov2CoytuXA2JYfrAwjMwXmkhPXtWhcxjZF4kugRIEHan+SuR4ytnWcfhnbGjEAzSKEEEKIxaIBXAghhMghPU2hl8oFpnd3OiWJjG3kwdg2+FgJRmeMalNJOHkytg2nsS1vzHmB7XOLMALL2NZxem1sa4QicCGEECKH9DQCL5eNXbu6HJ3J2DZgxrb8YGVYMQP1TSUytqUMm7EtL5S8wPTsEo3AMrZ1nN4Y2xZmAJpACCGEEItFA7gQQgiRQ3qaQqdslHaOsqsX+5KxbSCMbXkimNicZm0kY1uF3BvbckZIobf5k7cytnWc7hrbFkYRuBBCCJFDeh6BF3YWq6w1PY3GZWyj18a2PGElGJsp02obydhWTa6NbTlhzgtM72szAk+Rsa3jdN7YtjADcLhCCCGEWCwawIUQQogcYu7efKlO7czsPkKmbVvPdto9DkDH0QqPcPcDu7j9jhH1eRs6t4NGN48jN/oE9aEDSt/60J4O4ABmdq27P62nO+0COo7hZVjaRMcxnAxLe+g42kcpdCGEECKHaAAXQgghckg/BvAtfdhnN9BxDC/D0iY6juFkWNpDx9EmPX8GLoQQQoj2UQpdCCGEyCE9HcDN7Hgz+6WZ3Wxmb+/lvpeKmR1iZlea2Q1m9nMze2MsX29m3zazm+L/df2uayuYWdHMfmxml8XPjzSzq+M5udjM8vproG2TR32CNLqckEb7zyDps2cDuJkVgY8AJwBHAKea2RG92n8bzAHnuPsRwNHA62O93w5c7u6HA5fHz3ngjcCNyef3A+e5+2HADuDMvtSqz+RYnyCNLguk0YFhYPTZywj8GcDN7n6Lu+8DLgJO7uH+l4S7b3X3H8XpGcKJ20So+4VxsQuBU/pTw9Yxs83AS4CPx88GHAd8IS6Si+PoErnUJ0ijywhptM8Mmj57OYBvAm5PPt8Ry3KDmR0KPBm4Gtjo7lvjrLuBjX2q1mL4IPBWKr8GsAGYdp//hYLcnZMOknt9gjQ65Eij/Weg9CkTW4uY2Wrgi8DZ7v5gOs+DlX+g7fxmdhJwr7tf1++6iO4gjYpBJ88aHUR99vLnRO8EDkk+b45lA4+ZjRJE9xl3/1IsvsfMDnL3rWZ2EHBv/2rYEs8GXmZmJwLjwBrgQ8CUmY3EO8jcnJMukFt9gjS6TJBG+8vA6bOXEfg1wOHRsbcCeCVwaQ/3vyTiM47zgRvd/QPJrEuB0+P06cBXe123xeDu73D3ze5+KKHtr3D3VwNXAi+Piw38cXSRXOoTpNFlhDTaRwZSn+7esz/gROBXwK+B/9HLfbdR5+cQ0jo/Ba6PfycSnn1cDtwEfAdY3++6LuKYjgUui9OPAn4I3Ax8Hhjrd/362C6502estzS6TP6k0cH4GxR96k1sQgghRA6RiU0IIYTIIRrAhRBCiByiAVwIIYTIIRrAhRBCiByiAVwIIYTIIRrAhRBCiByiAVwIIYTIIRrAhRBCiBzy/wEq1lk2qj+Y/wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 504x216 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "n = 50  # nb samples\n",
    "xtot = np.zeros((n + 1, 2))\n",
    "xtot[:, 0] = np.cos(\n",
    "    (np.arange(n + 1) + 1.0) * 0.9 / (n + 2) * 2 * np.pi)\n",
    "xtot[:, 1] = np.sin(\n",
    "    (np.arange(n + 1) + 1.0) * 0.9 / (n + 2) * 2 * np.pi)\n",
    "\n",
    "xs = xtot[:n, :]\n",
    "xt = xtot[1:, :]\n",
    "\n",
    "a, b = ot.unif(n), ot.unif(n)  # uniform distribution on samples\n",
    "\n",
    "# loss matrix\n",
    "M1 = ot.dist(xs, xt, metric='euclidean')\n",
    "M1 /= M1.max()\n",
    "\n",
    "# loss matrix\n",
    "M2 = ot.dist(xs, xt, metric='sqeuclidean')\n",
    "M2 /= M2.max()\n",
    "\n",
    "# loss matrix\n",
    "Mp = np.sqrt(ot.dist(xs, xt, metric='euclidean'))\n",
    "Mp /= Mp.max()\n",
    "\n",
    "\n",
    "# Data\n",
    "pl.figure(4, figsize=(7, 3))\n",
    "pl.clf()\n",
    "pl.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')\n",
    "pl.plot(xt[:, 0], xt[:, 1], 'xr', label='Target samples')\n",
    "pl.axis('equal')\n",
    "pl.title('Source and traget distributions')\n",
    "\n",
    "\n",
    "# Cost matrices\n",
    "pl.figure(5, figsize=(7, 3))\n",
    "\n",
    "pl.subplot(1, 3, 1)\n",
    "pl.imshow(M1, interpolation='nearest')\n",
    "pl.title('Euclidean cost')\n",
    "\n",
    "pl.subplot(1, 3, 2)\n",
    "pl.imshow(M2, interpolation='nearest')\n",
    "pl.title('Squared Euclidean cost')\n",
    "\n",
    "pl.subplot(1, 3, 3)\n",
    "pl.imshow(Mp, interpolation='nearest')\n",
    "pl.title('Sqrt Euclidean cost')\n",
    "pl.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Dataset 2 : Plot  OT Matrices\n",
    "-----------------------------\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfAAAADQCAYAAAD4dzNkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO2de5wcVZn3f08IJMwkM4R0REDCRSLK7rIXgiv7siZjEAGDYMQxuC9JQHck++q6GgMEDTMJaExw0HXh5SKSAVkIg4ByCZfEalEUZcIKIhkDAY0J1yRcJzO5zTz7x6meqempvlR3Vdep7t/38zmfma6urq7Lt55T59TTp0RVQQghhJBkMSruFSCEEEJIcFiBE0IIIQmEFTghhBCSQFiBE0IIIQmEFTghhBCSQFiBE0IIIQmEFbhliEiHiFzu/v/PIrKhmHkJiQoRaRORWyr8nTwPSEUQERWRo93/rxWRxcXMawM1VYGLyDwReVpEekXkFRG5RkQOcN+7VkR63LJbRPZ4Xj/gs6zpIjLgmSdTTgxrfVX1l6p6TFjLI8EJ05lqhOdBsrDVZ9ejLQXm6XDXy+vZU2Guh6peoKqXhbnMKKmZClxEFgBYDmAhgEYAHwJwOIA1IrKfe+DGqeo4AN8CcHvmtaqelmOxL3nmyZTHKrJBJHIiciZWRGR0BIvleZAAbPU5oJMrsjz726jWKwnURAUuIg0AlgD4kqo+qKp7VPXPAJoBHAHg/0bwnX8WkZM9r4d1Q4rISSLyaxF5U0Q2i8g8n2UMuyoVkb8Xkf8RkXdE5HYAY7PmnykiT7rL/LWIHOd572IRed797HoR+aTnvXki8qiIfEdE3hCRP4mIlRVQpYjCGRFJich97vF5XUR+KSKj3PeGHVsRWeXpQp4nIo9mLcvb7fdxEfmdiLztutTmme8Id97PichfADju9A95/HtKRKZ7PnOkiDzirssaAKmg2+pZFs8DC7DM5+kiskVELhKRVwDcBuABAId4WtaHBFyXES14r3siso+IXOI59k+IyGE+yxl2O0ZEForIyyLykoicnzXvGNeVv4jIq2J6MPZ335vg7putrkv3ich7PJ/9uYhcJiK/ctfnYREJfJ7VRAUO4J9gTvK7vBNVtQfAagAfreTKiMjhMML+F4BJAP4OwJMFPrMfgJ8A+BGAAwHcAeBTnvf/HsCNAL4AYCKA6wDcIyJj3FmeB/DPMFfeSwDcIiIHe77iHwFsgAnWKwD8UESkrA1NNlE4swDAFphjfhCASwBooWNbBDsAzAFwAICPA5gvImdlzTMNwAcAfExEDgVwP4DL3e/7GoA7RWSSO++tAJ6AceEyAHODbWZx8DyoKLb5/G73vcNh3D0Nw3tyXiphffLxVQDnADgdQAOA8wH05vuAiJwKc258FMAUACdnzfJtAO+D8fZoAIcCuNR9bxSAlTDbNxlAH4Crsj7/WQDnAXgXgP3c7wpErVTgKQDbVHWvz3svo/QWxiHu1ae31Bfxuc8CWKuqt7lXwttVNW/ggunu2hfA99zP/BhAl+f9FgDXqepvVbVfVW8CsMv9HFT1DlV9SVUHVPV2AM8B+KDn85tU9Qeq2g/gJgAHw5yUtUoUzuyB2a+Hu8fwl2oeRlDo2OZFVX+uqk+7x/b3MC2aaVmztanqDlXtg2ltrVbV1e5n1gBYB+B0EZkM4AQAi1V1l6r+AsC9BVaB54H92ObzAIBW17G+AN/5tSzPbiryc58H8A1V3aCGp1R1e4HPNANYqap/UNUdANoyb7gXdS0AvqKqr6vqOzC3HWYDgOvynara6773TYw8J1eq6rPu9nfCXAgEolYq8G0AUuJ/r+Vg9/1SeElVD8gqO4r43GEwLYEgHALgRfcEybDJ8//hABZ45Xa/5xAAEJE5nm7FNwH8NYaftK9k/lHVzJXpuIDrWE1E4cwVADYCeFhEXhCRi93phY5tXkTkH0Uk7XbXvQXgAowMyJs9/x8O4NNZrpwEs12HAHgjy+NC68LzwH5s83mrqu4s4Tu/k+VZsb1DpbrmPW+82zAJQB2AJzwuPehOh4jUich1IrJJRN4G8AsAB4jIPp5lvOL5vxcleFYrFfhjMFfhs7wTRWQcTNfNzyL4zh0wBzjDuz3/bwbw3oDLexnAoVndeZOzlvnNLLnrVPU2t6vyBwC+CGCiqh4A4A8AqqFrMCpCd0ZV31HVBap6FIBPAPiqiMxA4WM7zCUR8boEmC7vewAcpqqNAK7FyGPrDaabAfwoy5V6Vf22uy4TslrQk1E6PA/swCafgeE++r0OSvY5sg/cytSlVNe898m927ANplv8rzyeNapJAATM7YVjAPyjqjYA+HBm1QKuQ15qogJX1bdg7nf9l4icKiL7isgRMN0WW2Du1YTNkwBmu981FcDZnvf+G8DJItIsIqNFZKKIFOo+eQzAXgD/7i5zFoZ3/f0AwAVua0xEpF5MctN4APUwJ8hWABCR82BaHiQHUTgjJrnqaDewvQWgH6YrsdCxfQrAX4nI34nIWHi68lzGA3hdVXeKyAdhuqbzcQuAM0TkY2KSe8aKSQJ6j6pugulOXyIi+4nISQDOCLqtHngeWIBlPvvxKoCJItIYdD1cngUw1j3W+wL4BoAxnvdvAHCZiExxvThORCYWWGYngHkicqyI1AFozbyhqgMwrn1XRN4FACJyqIh8zJ1lPEwF/6aIHOj9bJjURAUOAKq6AibJ4jsA3gbwW5irshmquqvExXqzJjMlk6yxGOaK7w2YE+dWz7r8BSaZYgGA12GCXN6fQ6jqbpir53nuZz4DT0KKqq4D8K8wiRJvwHRtzXPfWw+gHebEehXA3wD4VYnbXDNE4MwUAGsB9MAci/+vqukiju2zAJa6n30OwKPDF4t/A7BURN6BSaLpLLBdmwGc6W7bVnebFmIoHnwWJpnrdZjAc3OB7eJ5kABs8TnHuv0RJnfjBbdLOlcW+oVZnm1zP/8WzHlwA4AXYVrk3qz0K2HOi4dhtv2HAPYvsE4PAPgezC83Nrp/vVzkTv+N202+FqbVDfdz+8O01H8D070eOjL8NgUhxAZEpAPAFlX9RtzrQki50OdoqJkWOCGEEFJNsAInhBBCEkgoXegiciOAmQBeU9URSSFiRnn6KYA/uZPuUtWlZX8xIQWgm8RW6CYpl7DGRe6ASRrJl+zyS1WdGdL3EVIsHaCbxE46QDdJGYRSgavqL9yfJIRKKpXSI44IfbEkJp544oltqjqp8JzhQTdJMdBNYjO5/IziyUS5OFHMo99eAvA1VX3GbyYRaYEZog6TJ0/GunXrKriKJEpEpOjRxSoM3axx6CaxmVx+ViqJ7X9gxsv9W5gHF/wk14yqer2qTlXVqZMmVfSCmNQmdJPYCt0kealIBa6qb6t56g1UdTWAfaWER6cREjZ0k9gK3SSFqEgFLiLvzoyL6w71OApAoSfBEBI5dJPYCt0khQjlHriI3AZgOszTbrbADL+4LwCo6rUw4x/PF5G9MOPDzs56Ug0hkUA3ia3QTVIuYWWhn1Pg/asw8mHmhEQO3SS2QjdJuXAkNkIIISSBsAInhBBCEggrcEIIISSBsAInhBBCEggrcEIIISSBsAInhBBCEggrcEIIISSBsAInhBBCEggrcEIIISSBsAInhBBCEggrcEIIISSBsAInhBBCEggrcEIIISSBsAInhBBCEggrcEIIISSBsAInhBBCEggrcEIIISSBsAInhBBCEkgoFbiI3Cgir4nIH3K8LyLyfRHZKCK/F5F/CON7E8eKFUA6DQBoa3OnpdPA6af7T1+xotJrWHXQzSKhmxWHbhYJ3cxJWC3wDgCn5nn/NABT3NIC4JqQvtdOcgn3/PNAczOQTmPJEndaczNw8skjp59xBjB69Mhl1JCcIdEBujlEGG42N5v5azx4hkAH6OZw/PwcPRqYOZNu+qGqoRQARwD4Q473rgNwjuf1BgAHF1rm8ccfr1azfLmq46iqamurO81xVFtaVFMpVcdRwJ3mvs78vwSLh6apjpze3u6/jJYW/+9cvrxy210iANZpSL4FKTXppqq/n+3tqnV1g169ebejeyekdMO1jv7uSkd7xxkH+8an9GffcPTWW1XXft3RHfVm+q7GlD5zlaMv3Ojo3gNT+tJ/e/ysrzfLV7pZbKGbRcZONx4WjJueGDtsGQl2UzW3n5US8T4AJ3le/wzA1BzztgBYB2Dd5MmTI94tZeInSp6KurXV7PElWKwK6BIsVkB12jT/6SvnFCmnV2iLsTRIVqebqiNceeMuR3cfkNKnz2sfrJB76lLaMdfJ6eCUKbmd7ZjraE/dUIX/3Px23TMhpTvuo5vFFrpZOHaunOMEiputrToy/uZqECXATdXcflonordYdSVZ4IqxmIraV6wcV5K5pPVdRkJa5UkPkt6SBDd3X75c11/tDKusv3ycv1cXXaS6d01hNwdSKe1b7eiCBf6Of/m4oUp9V2NKt//Y0YFv5zh3LPKTbkZErriZmV5k7PRtzGSWVc4y8q2fRcRdgSe/K8jnirGnLqXT4RRfUee66sx1dejTbeQn53SYdbH96tLSIFl1br5wo6N941P6lb/zr6xvPDcEN31aSjf8S46Lgw86uqsxpZs67PWTbkZEDq8CNVAKxMdyWvGZSj2psbNSIn4cwAMABMCHADxezDJjEzFAaztQRZ2rpXzaaQXvVRa8j+43zbKrS0uDZFW4uXNuy+D96566lN5+gaMPPKD62u2ODpSTW+HnpudcyHfROTAxpX+8xtGbblK9ad5Qy3zPhJT2/6tdvUZ0MwSCxM3Me+Xcvy7DzbzrYWGPZqQVOIDbALwMYA+ALQA+B+ACABe47wuAqwE8D+DpYrqBNE4Ri2xt52z5hiFAgCSPXD0Btl1dxhEkq93N7T929K39Ahz/sJJ5ikiQy9famg5H3x6b0j0P081qdbOk2BRV7Azopm09mpG3wKMosXYFFdPKjeNKLYyr3JiIq5UTRYnbzd0HDLW2b/mco48vd7R/Ysw9MEXe79wzIaX3LXCGJcEN0M2qcdPK3sGA9+Jtipuquf2MXbZ8pSIi+hzYlXMcXYjl1l+VDZKjy97v6nIhlpurX6189xCDZEB83PzBZxPmpmrRLTO6mWw3Exc3VRPhpmpuP2OXLV+piIg5Kr8R90ssvC8ySJCryxh/TsEgGZAsN5+73mSUP/ixdt3ZkBA3VYvqNeqtT+mDp7Rr3/iUvnMP3SynMG4GoNgezZh/hsYKPB8WVXKhUuxJFrOESSyVdHPPhKHu8sc+0679B1afm3vXOLq7MaUPn9YeS7c63SwBxs2KrRIr8Bz4/Swr7u6S0Ciim2sNZuh0mJ+nRb2tDJLBoJt0s5RCN8ukyNsDldxeVuCqwRPAqhXvlXNjo2pDQ0WumhkkC+Dxc/Fi1fvvV33wlHbdvV9dTbmZ+dlb35hG7RvboK/dTjeDlKjdzM7sriU34+pxYAWu6t81UsEKzAry7IOoT0QGyQK4x+atn5hj8+Ap7TogonuvaK9JN1+73dGdYxu1d0xD5F3qdLMA2XGjvV1VRLW9Nt2sdLc6K/AM2VdRtidZhE3WlXRrq+p0OLoGM9T7O83BfREiDJKFefa6oZ9W9e+f7AcwBIZuWu3msNiZ8IeDBKbIbvVK+xm7bPlK2CLmHaO8lqnQickgmZ//+A/6OQJPl3pPXUpf+GiL9q8N/4KbbuaHsTMHMcfO2GXLV8oWkfdtClPBrjEGSQ9Zbj7/vOotn3P06X9qGTkgS62S5eaT33W0b0yj7tq/QfeuCfeWF93MgrGzMBbEzthly1fKFrHW79sUQwVPVAZJD1lu3nyeozv3b9SBWsrHKISPm+uvdrTr+KFx33lxGYGbqoydxWBB7IxdtnwlTBFr8r5NCUTZVcYgmYVjntK1BIu1d1xKd59XY/kYJRCVn3TTB8bOQMQRO2OXLV8pV0TetymR7EQ/tsDppk14Lnx2NbAFnl3CqMDpZ4lUOHbGLlu+ElhE3rcpn+yus5YW8zMzxyn7yrvmg6THz0suUb3zi472jW0wWdZ0szg8fi7EcvM7+QNSgwNq0M3y3WTsLJEYYmfssuUrgUXkfZvyyT6RHcf8Tr6lpex9WPNB0t13e9c4Oh2O9o5p0P7xjTodDt0sFo+fN57raN94Mzb8QiynmyG4ydhZBjHEzthly1fKEZH3bUIkpG6hmg+Sqtq/1tHecSldgxm6d1w4V+e1zO6HHO2tH3pMKd0sowudsTN8Io6dscuWrwQVkfdtwifMfVrrQZJ+hgvdpJs2Uwk/Y5ctXym7Bc4un3BgCzwUN3t6VG9rcUL/CVRN43Gzpy6lT32PbobSAqeb4RBx7ByFEBCRU0Vkg4hsFJGLfd6fJyJbReRJt3w+jO8FAKxYAaTTAICOuWmguRlYtAg9GAd0dprX7vukBNLuPu3sNPt00SKgudns68z7K1bEu44FiM1Pj5ttbUDXijTOvOks4JOz0Iql9LNcPG62YileOXEW3nfRJ/H6nWm0tXnmsdhPxs4qphKx069WD1IA7APgeQBHAdgPwFMAjs2aZx6Aq4Iuu6gryays1Mzg8uVmpRIXT2LGyjlDA/iXkjSEGFo5UfkZ1E1AtesfWnRPPe97h0ZW0lDfapPV//Q/BU8aqjk3VRk7o6YCsTMMCU8E8JDn9SIAi7Lmia4CV2XXTyUpY1/HFCQj8TOIm5mhUXvHpQbH8SbR8Mdrhx4Gk4CLy3jdVGXsrCQRxM4wutAPBbDZ83qLOy2bT4nI70XkxyJyWAjfC8B0TcpHmrB023xcisuwdNt8yEeahrrQSGgkdF/H5mdmf12+3eyvK3rmY5+Trd9fiaWtDXj/BU1o702Mn4ydNUJk+9qvVg9SAJwN4AbP63ORdcUIYCKAMe7/XwDg5FleC4B1ANZNnjw58isbEpDktcBD87MUN3vuLa1FSErE0+Oxq9H6Fnisbmb2F2NnhbC0Bf4iAO9V4XvcaYOo6nZV3eW+vAHA8bkWpqrXq+pUVZ06adKk3N+aScBIp7FjpkkUSKMJmDWLyRdRkcyEttD8LMXN/rObcWcz3awIrp+jPjULaTSh81Od6Pl48+CxoJtg4lpcRBQ7w6jAuwBMEZEjRWQ/ALMB3OOdQUQO9rz8BIDusr/1hBPMDlm1CjN7OwEA99U1A7NnGxG7usr+CpJFV5fZt01NOHbOCcCyZcCiRVh/c9eQoCecEPdaZlN5P103d920Cmfu6sSUKXSzImT8nD0b9+7fDFXgjL5OYNUqupkhEzfTaXPeLloELFtmzuemJvoZFVHFTr9medAC4HQAz8JkVH7dnbYUwCfc/5cBeAYmyzIN4P3FLLdgMga7f+Il4P5HTL+1jcLPYtzM/N67fyLdrDiOGWZ1CRbr3iJGaKs1Nxk3Yyak2FlxYYOUfCJy5KB4KWX/xxUkoyh0017oJt20mTD9jF22fIVXkpaTkBZ4FKWQm09/n8lrseJxc0e9eYBMPmrJTcZNC6j1FviIp+fwaTmVJbO/W1oGn6bVU+fu/xwDQNRKkBz4mek+v38h3YyFrAFKHjylXfdMyD9ASa24ybhpASHGzlCGUq0YngzKNd9ykwIAPHzyCiZgVBpPwtB9dc0AYJIJ7U0Yih7Xz7fWdmHRUZ04+JwmrJyTHtpXdLMyeBKGPnDuCfjnXy3Dk6dZn2wZLZnY2dWFjtPNvqGbMRFm7PSr1W0pI64kefVoJ0V2B6HaWznufkhfavzc/RD9tIHHvlX4dkatuMnYaRllxs7YZctX8onI+zd2ECQho+qDpKr23s973zZRrJ+14CZjp12EETtjly1fyRaRGZSWwha4qtJPa3EKP8KVbpJYYAucV5GxEqBrrtqDpKrqE99hC9wqXB9/+23j556H/f2sBTcZOy0jhNiZnCS2FSuAK680w8+dbp7/i0WLgJkzOfRfnHiSYFbOSQNNTSZJpqvL1uErw8czdOp7L2nGL/4fh061hq4uYNYsHHooMG0a8NV73WTXVatqw02AsdNWXDcBoLUVaHukBDf9anVbyrArScdRra9XbW9XYPjzVfnMWgso4moS1drKcbd1z/nmZyG/u7Lwz0JIBXEcHUil9Efn16Cb7vYzdlpKka3wXH7GLlu+kisLnV1AllLg+FR7kNzV6A7deSDdtA7H0R31Neqmu/2MnZZSxLFJfAXOJAy7Keb4VGuQpJt2U8tuFrv9JB7K/ZVE7LLlK2yBJ4wab4FnHp5BNy3EcXTPhNwPl6lqN93tZ+y0lFpogXMgAsvJOj6D99lq4T6ju+13fYluWol7PLbdYY7Pd8+sITc928/YaSFFxE3V3H4mIwt9xQqTmdfZibZHmtDa6k6fNYtDANqCZ/jK1lbgvJtrZGhb183+2zrx/aebcO657nS6aQ+umw1nNmHaNOArP60RNwHGTtspN2761eq2lMErSV5BJorf/EbNccoC1djKcV18827j5nM/oJs285//WUNuqjJ2JohXXzVu7to18r1cfiajBZ55UElzM5bgUvPbWveqhdhDWxsgAnzoQ+a1iCltbXGuVcS4bu4/z7h5xIV000Yybn75y+Z1TbgJMHYmgIybBx1kXo8ZE8BNv1rdlpK5kmQWZXJ4/HHVtrbaaeXQzWTR01M7bqrSz6QwMKD6ve/5u6ma289EtMDb2gB10rg0dQ2WYjEuTV0DddLVf/WcQNavB1KpuNeicmTcXHSAcXMx3bSa+vq416CyMHYmg1deAd58M/jnQqnAReRUEdkgIhtF5GKf98eIyO3u+78VkSMCfUHmGb6d7jCAbpcQhwG0ix07gE2bgA98AEPJMhYQqZ+um91txs2B2+im7dSMmwBjZ0JYv950m19ySbDPlV2Bi8g+AK4GcBqAYwGcIyLHZs32OQBvqOrRAL4LYHmgL7niCjN2b5PJIkVTk3l9xRXlrj4JkQ0bAFXg2GPtubcYuZ+um9v+pgmHHw6MmkE3badm3GQGeiJQBbq7gSOOAL75zWCfDaMF/kEAG1X1BVXdDWAVgDOz5jkTwE3u/z8GMENEpOhvWLgQWLYMSKfxyCMwV4/LlpnpxBq6u4EJE4aSMSwhWj9dN+sfT2PTJkB+TjdJ0UTr5gknAHfdBQBYsgRom+a2xmfPBi68MJQNIOWzdSuwfbvpuQxKGBX4oQA2e15vcaf5zqOqewG8BWCi38JEpEVE1onIuq1bt5qJzKS0np07gRdeMBIGuDSrBKH5mc/N4y6nmyQwFXGTcdNuurvN37gq8FBR1etVdaqqTp00aRIAN83+I01Yum0+LsVlWLptPuQjTdZ0hRHTfT4wUJqESSGfm1fsoJskPhg3k0t3NzB5MjBuXPDPhlGBvwjgMM/r97jTfOcRkdEAGgFsL/YLmElpP93dQEMDcGh2+yF+IvUz4+bFjcxCJ4GpiJuMm/ayfTvw6qulN3zCqMC7AEwRkSNFZD8AswHckzXPPQDmuv+fDcBxf9tWHF/4AnDWWcMzKc86y0wnsbNrF7Bxo5Xd50DUfqbTwMyZeHneIrRiKX54ipvle+WVJomIkNxUxE0sMm52nE43baOc7nMghArcvS/zRQAPAegG0Kmqz4jIUhH5hDvbDwFMFJGNAL4KYMTPJQqSXTNYWFPUKs89B/T329l9HrmfXV3AZZfh0JuX4cvHpfGvt7pZ6JdeapKICMlBpdzEsmVYOSdtxtmmm1bR3W16LRsbS1yA3+guthS/MX35ODz76OxUveIK1f7+/POhSke7UlV94UZHe+roZ1KpZjcZO+3kjTfMqJWPPlp43lx+WpfE5geTMexlzx7TAn//+4FRibApfNragKPOb0J7L/0kdsHYaS/ldp8DFmah+8FkDHvZuNFU4sdmDz9RQ7S1Ae/ck8aCOvpJ7IKx0166u82YGQceWPoyElGBczhAe+nuBvbfHzj88LjXJEbSaYw7vxmr59JPYhmMnVbyzjvA5s3lN3ySUYF3dZnh/2DGMW57xB2gYNUqZlPGyN69wLPPAsccA+yzT9xrEyOunxMmAKefTj+JRTB2WkkY3edAUirwCy80w/81N6NtWhpLlrjT77qL2ZQx8qc/mZ+Q1XL3OYBBPz98dTM+8y76SSyCsdNKurvNUxvdMXdKZnQ4q1MBhg0LOB9ovobDAsbM+vXm4fNHHhn3mlhAUxOeu7wTn7qoGX+mn8QmGDutIvPUxpNOKn9ZyWiBg9mUttHfb4ZPfd/7gNHJuQyMjLY24K+/xEx0Yh+MnXbxxz8OPbWxXBJVgTOb0h42bQL6+uwcvCUO2tqA/rVpLKynn8QuGDvtIsynNiamAmc2pV2sXw/suy9w9NFxr4klpNMYNbsZT33d+Ln3VvpJLIGx0xr6+kzuUFjDTienAu/qMuJ1dWHlnDTQ1GTG9u3qMiIyo7JiDAyYbqApU0wlTjCY7XvYYcC0acBX7mG2L7EEZqJbw7PPmvgZVuJvcirwCy80SRcnnIB5q83V43k3m9dobmZGZQXZvNkkYrD73IOb7XvIV5rx+femcdVV7nRm+5K4YSa6NWSe2njIIeEsL3npR8yojJ3ubvO77ylT4l4Ty2hqgnR2YtYZzXgB86HN10DoJrEBxs3YyTy1cerU8J7FlZwWuAszKuNF1VTgRx9tfkJGhsi4+Z0dxs3L6CaxBMbN+IniqY2JrMCZURkfL74IvP02u8/9yLi52HXzwvF0k9gB42b8dHcD48YBhx0W3jITV4EPZlTOmoU0TLfQjpluRiWT2SKnu9s8deyYY+JeEwtx3RTXzdtndaLn43STWIAnE70H48xzwZub0TE3PfQ+/YyMqJ7amLwKPJONPns27qtrBgDM7HUzKpnMFimZ7vOjjgLGjo17bSzE4+a9+zdDBDijj24SC8i42dSEY+ecACxbBixahPU3dw1V7vQzMjJPbQy75zJ5SWwXXjj4b/19JimjCfNNRiWTMiLl1VeBN94IZwjAqsTj5rj7O/HpT5hhVfWuu5jMRuLF4+a8m5qAtImd4zB/sGVOP6Mj89TGI44Id7lltcBF5EARWSMiz7l/J+SYr19EnnTLPeV8ZwYmZVSe9etN9mRSus/j8jPj5hU9TGYj/jB21g7epzaG2X0OAFDVkilAKQsAAA4YSURBVAuAFQAudv+/GMDyHPP1lLL8448/XvPiOKqplC7BYtVUyrwmkXHVVaodHaV/HsA6LcO3oCVKP4txc8B1c0d9SvtW002bqSk3VRk7K8iGDaptbarPPlv6MnL5We71wJkAbnL/vwnAWWUur3iYlFFRtm4Ftm1L3KND4/Ezk8zmuvnLkxZh1Gy6SYYRf+xkInBF6O6O7qmN5VbgB6nqy+7/rwDINTz7WBFZJyK/EZFwRGVSRkVZv978ff/7412PgMTjZ5ab0x9bhvSJdJMMI/7YyUTgyMk8tfGYY6J5amPBRYrIWgDv9nnr694XqqoiojkWc7iqvigiRwFwRORpVX0+x/e1AGgBgMmTJ+deMSZlVJTubvP7xfHj416T4VTSz1Ld7L2/Eyc1N+O3dLOmsNJNgInAFSTypzb69asXWwBsAHCw+//BADYU8ZkOAGcXs/yi7uWoamurKqDmfg6gS7BYATOdlM/27eYezq9/Xd5yUPn7jJH5STerC7pJN6Pg3ntVv/lN1d27y1tOLj/L7UK/B8Bc9/+5AH6aPYOITBCRMe7/KQD/B8D6Mr93GBxlKFoy3ecJHH0tdj+zR2f7Wv012HEf3ST2uMm4GQ2VeGpjuRX4twF8VESeA3Cy+xoiMlVEbnDn+QCAdSLyFIA0gG+raqgVOBPaoqW72zw954AD4l6TwMTvZ1ZC268+vAijzmk2j8TNvE83axFr3Bx8TvisWcBZZwFpTyVOP0umEk9tLOu2uqpuBzDDZ/o6AJ93//81gL8p53sKMixpCMCyZpPQtqALmIeh+44kMG+9Bbz0EjBjxFG2Hyv8zHJz+j3NcE5chO4fdQHngW7WKLa52doKYNps4PbbgVWrsOT6JrRNS9PPMli/3iSuRfnUxuSNxOYHE9oio7vb/E3Yz8fsIctNdTrx4TNNQtveTzVj9J10k8SEx03T4m4C7r7bfeToQXzkaBmomu7z97432qc2Jm8s9AJwlKEQWLHCdJ3BVOAHHQQc+BS70sqlrQ0YNWNohLZvvUE3S8Lj5yDs6i0bxs4QcN0c9tTGKN30y2yzpRSbTTkC7yhDdXWq7e2q6smudBzV5ctLW3Yt4O6/Hfc52tam+uR3nVBGa0KFM32jLGG4uWvfOv3lrHbt66ObgXD3Ye/9zrDX5fhJN10YO8vD3X+PL3d06VLVXQ9GGztjly1fKUlEz8kMqBFQRLW93bwO4WSvCRxHdx+Q0p9/eLH2Twxnf9V8kMxyc9sl7ToA0cdn082gvHiLozvqUrp1fjhDgda8m6qMnSEx8DNHd9Sn9PdnhTdMbe1U4MuXD+6wwavG9nbV+nqO+xuQN//d/D504BuLQ1lezQdJHze3LGjXXfsaNwfoZlFs3Kh6+eWqXacbP3Vx+X7WvJuqjJ0hsXev6l/mheemam4/Y5ctXymrK8iFgxUUQZ4TVxdHfxWZxEI3K0SWm889p3rzeY4+c1KLueAJyU+6ORL6WQQxx87YZctXwhKR93UKkKfrLPv9cmCQ9MHz1LJd+9bp4+e0686ddHOQLDdvPs/RnWMbdKChcchH3gOPxk1Vxs5CxBw7Y5ctXwlFRN7XKY4cJ+qw98s8URkks8hyc8sCc0/8V2fTzWE4ju5qNG72jkvpnvNbRu6TMv2kmz4wdhZHjLEzdtnylVBE5H2dkWTtk9ZW1elwdA1maJRdZQySWfi4+erFQ/fE9x5INxcvpptBS2gVOGPncHz2x8o5ji7E8shvM9RuBZ4F7+voyCtrx1FtaFBtbIz0xGSQzA/d1BFurl7oaO+YBt1dTzeLLVG4qUo/feNmKmUuajIt8Ar7Gbts+UpUIg7r8kilVFtaRl5pVvu9He8+cCvvEWKGLCKDZBF4jkvfuJR2Hd+iv1ji1NZ9cfcnjEuwWHvHNOjecXQzSInMTVXeE8+uO9zKOy4/Y5ctX4lExDytz6gPQiwU0e2zBjN0OpzB7nRVjeREZJAsQJab/Wsd3V3foH1jGvWOf6sNN2/4F7pZbon64rIm7okX2V2+EMt15ZzoG3+swDP43ddxHNMKr8ZWeYzdPtkwSBYgh5tvzW7RnjpzrHY2DE/iqiY3N3WYATAePKVddzbQzVJLZBV4MffEGTcjgRV4Hvzu7UyHoz110XeNVISYun2yYZAMTi433x6b0v611eHm3gONmz11KX10Vvtg8h7dTKabjJvhwwq8ENkHy3H8p9mMZd0+2TBIlojHw/4DU/rgRY52zHW0tz6VnBHcfNy8bjbdjKIwbgbE8ripmtvP2GXLVyomok93SU9dSqfDseYAFoVl3T7ZMEiWQAA3bzw3OW4+epnZjgc/1q594+lmmMU2NwHVh0/OcevSBj8tj5uquf2MXbZ8pWIiFntfPFcXSqXv++Ra38x0S7p9smGQLIE8bmZGcOsdZ+4b76hP6f+02+lm/7Llg/e4M93l6z/fPvigHLpZXW7mapUPO9ZxZLGXG+stu8CMXbZ8paJdQdkEuSrLNW8YBzvfWLtZ37dyjjPinpRNvQYMkiGR5dvAzxzdOyGl6z7bPpjs1jsupRuudcwjN2N087XbHd3VkNIF/0A3K1VscnOYb5VsXBSoqJPW2xpJBQ7g0wCeATAAYGqe+U4FsAHARgAXF7v8WEUMOupOtpy5Wj6nnRZsuo9webt3LL7/VOkgGaWfSXFz+nTVBy40leiw++V+Aa6lxRQt3c3+iSl9aeFQt3hPXUpvPs/Rzk6TZT5AN2vOTXUcffjk5b6DwKyc4xOvclW+fh7mcjZX3PS7kMg1zSKiqsA/AOAYAD/PJSGAfQA8D+AoAPsBeArAscUsP1YRc+FzoANlY+a66sx3NZr1nX4tbe/JYEu3TzYxBMnI/LTdzYFUSrfd4eicOSPdBFQvn+HongkpffGWIsZD8HFzIJXSN+92dPPNzuBPvnrrU/rl4+hmMaXm3FQtP475xch8Y3iUG6stcVM1ogp8cCH5JTwRwEOe14sALCpmudaJGKR7KN9VXYDpuYYvDHTlakOiiFY+SGZKFH4mzc09E4Z+R37/1xy9/HLVjrnOYJf7jrqU3vlFR+9f6GjfuKF5H2lz9I47VJ3Fjva603fUp7RjrqPTpvm7efXZPi1tulm7bqoGuyXpmb+o6QHiZqDeUkvcVM3tZyUkPBvADZ7X5wK4Ks+yWgCsA7Bu8uTJEe+WgJR5XwXQnEEv13Rf4SxLsCgWS4Nk0X5Wi5uaSumN5/q3fvymHX20//QvfEH15VuHfsdNN+lmTsrM5ckXI4uOm57Kvlr8LEawtQD+4FPO9MwTWgXuLVZeSfoRJBsz817Aq8thwiV0DOIogmRcfibeTb9fLATtNaKbdLMcgrqZeb8YZ3NV1AlobfsRZwu8errQg1Cgy6jo6QkVzg9LWznV000ZBD8/A9wDp5t0MzKCxM5czlaRm6q5/ayEhKMBvADgSAwlYvxVMctNtIhBMinzTU+ocH5YGiRL8jPRbqqGk4VON+lmFISRhV5FbqpGVIED+CSALQB2AXg1c7UI4BAAqz3znQ7gWZiMyq8Xu/zEi0iGUekgGaWfdLO6oJvEZnL5ORploKp3A7jbZ/pLrniZ16sBrC7nuwgJCv0ktkI3SRiMinsFCCGEEBIcVuCEEEJIAmEFTgghhCQQVuCEEEJIAmEFTgghhCQQVuCEEEJIAmEFTgghhCQQVuCEEEJIAmEFTgghhCQQVuCEEEJIAmEFTgghhCQQVuCEEEJIAmEFTgghhCQQVuCEEEJIAmEFTgghhCQQVuCEEEJIAmEFTgghhCQQVuCEEEJIAimrAheRT4vIMyIyICJT88z3ZxF5WkSeFJF15XwnIcVCP4mt0E0SBqPL/PwfAMwCcF0R8zap6rYyv4+QINBPYit0k5RNWRW4qnYDgIiEszaEhAj9JLZCN0kYVOoeuAJ4WESeEJGWfDOKSIuIrBORdVu3bq3Q6pEapyg/6SaJAbpJclKwBS4iawG82+etr6vqT4v8npNU9UUReReANSLyR1X9hd+Mqno9gOvd794qIpuyZkkBqLbupFrZpsPD/pJK+kk3qwa6WT1U43YV7WfBClxVTy53bVT1RffvayJyN4APAvCtwLM+Nyl7moisU9WcSR9JhNtUOnH5STeTC92sHqpxu4JsU+Rd6CJSLyLjM/8DOAUmgYOQ2KGfxFboJilEuT8j+6SIbAFwIoD7ReQhd/ohIrLane0gAI+KyFMAHgdwv6o+WM73ElIM9JPYCt0kYSCqGvc6BEJEWtz7PVUDt6k6qMZt5jZVB9W6zdW4XUG2KXEVOCGEEEI4lCohhBCSSFiBE0IIIQkkkRV4seMI246InCoiG0Rko4hcHPf6hIGI3Cgir4lITWbLVoubQPX5STfpps2U4mciK3AMjSNc8LfktiIi+wC4GsBpAI4FcI6IHBvvWoVCB4BT416JGEm8m0DV+tkBukk37aUDAf1MZAWuqt2quiHu9SiTDwLYqKovqOpuAKsAnBnzOpWNO0rU63GvR1xUiZtAFfpJN+mmzZTiZyIr8CrhUACbPa+3uNMIsQH6SWyFbrqU+zjRyAhpHGFCQoduEluhm7WFtRV4GOMIW86LAA7zvH6PO41YTg24CdDPREI3awt2ocdHF4ApInKkiOwHYDaAe2JeJ0Iy0E9iK3TTJZEVeK5xhJOEqu4F8EUADwHoBtCpqs/Eu1blIyK3AXgMwDEiskVEPhf3OlWSanATqE4/6SbdtJlS/ORQqoQQQkgCSWQLnBBCCKl1WIETQgghCYQVOCGEEJJAWIETQgghCYQVOCGEEJJAWIETQgghCYQVOCGEEJJA/hdoiNdF+hXfnAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 504x216 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "G1 = ot.emd(a, b, M1)\n",
    "G2 = ot.emd(a, b, M2)\n",
    "Gp = ot.emd(a, b, Mp)\n",
    "\n",
    "# OT matrices\n",
    "pl.figure(6, figsize=(7, 3))\n",
    "\n",
    "pl.subplot(1, 3, 1)\n",
    "ot.plot.plot2D_samples_mat(xs, xt, G1, c=[.5, .5, 1])\n",
    "pl.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')\n",
    "pl.plot(xt[:, 0], xt[:, 1], 'xr', label='Target samples')\n",
    "pl.axis('equal')\n",
    "# pl.legend(loc=0)\n",
    "pl.title('OT Euclidean')\n",
    "\n",
    "pl.subplot(1, 3, 2)\n",
    "ot.plot.plot2D_samples_mat(xs, xt, G2, c=[.5, .5, 1])\n",
    "pl.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')\n",
    "pl.plot(xt[:, 0], xt[:, 1], 'xr', label='Target samples')\n",
    "pl.axis('equal')\n",
    "# pl.legend(loc=0)\n",
    "pl.title('OT squared Euclidean')\n",
    "\n",
    "pl.subplot(1, 3, 3)\n",
    "ot.plot.plot2D_samples_mat(xs, xt, Gp, c=[.5, .5, 1])\n",
    "pl.plot(xs[:, 0], xs[:, 1], '+b', label='Source samples')\n",
    "pl.plot(xt[:, 0], xt[:, 1], 'xr', label='Target samples')\n",
    "pl.axis('equal')\n",
    "# pl.legend(loc=0)\n",
    "pl.title('OT sqrt Euclidean')\n",
    "pl.tight_layout()\n",
    "\n",
    "pl.show()"
   ]
  }
 ],
 "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.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}