From 631e96e41a8e9a546503577d2a484d79748ee05b Mon Sep 17 00:00:00 2001 From: ritchie46 Date: Sat, 5 Oct 2019 18:50:17 +0200 Subject: [PATCH] working 2d change of variables --- .../change_of_variables.ipynb | 194 +++++++++++++++--- 1 file changed, 166 insertions(+), 28 deletions(-) diff --git a/bayesian/normalizing_flows/change_of_variables.ipynb b/bayesian/normalizing_flows/change_of_variables.ipynb index e0dff56..5c1051b 100644 --- a/bayesian/normalizing_flows/change_of_variables.ipynb +++ b/bayesian/normalizing_flows/change_of_variables.ipynb @@ -2,7 +2,7 @@ "cells": [ { "cell_type": "code", - "execution_count": 3, + "execution_count": 1, "metadata": {}, "outputs": [], "source": [ @@ -28,7 +28,7 @@ }, { "cell_type": "code", - "execution_count": 221, + "execution_count": 2, "metadata": {}, "outputs": [ { @@ -37,7 +37,7 @@ "1.0" ] }, - "execution_count": 221, + "execution_count": 2, "metadata": {}, "output_type": "execute_result" }, @@ -75,7 +75,7 @@ }, { "cell_type": "code", - "execution_count": 222, + "execution_count": 3, "metadata": {}, "outputs": [ { @@ -84,7 +84,7 @@ "5.641895835477563" ] }, - "execution_count": 222, + "execution_count": 3, "metadata": {}, "output_type": "execute_result" }, @@ -137,15 +137,14 @@ "\\end{eqnarray} $$\n", "\n", "$$\\begin{eqnarray}\n", - " P(y) &=& P(x)\\frac{\\text{d}f^{-1}(y)}{\\text{d}y} \\\\\n", - " &=& P(x)f'^{-1}(y) \\\\\n", + " P(y) &=& P(f^{-1}(y))\\frac{\\text{d}f^{-1}(y)}{\\text{d}y} \\\\\n", " &=& P(f^{-1}(y))f'^{-1}(y)\n", "\\end{eqnarray} $$\n" ] }, { "cell_type": "code", - "execution_count": 223, + "execution_count": 4, "metadata": {}, "outputs": [ { @@ -154,7 +153,7 @@ "0.9987379589284238" ] }, - "execution_count": 223, + "execution_count": 4, "metadata": {}, "output_type": "execute_result" }, @@ -194,16 +193,16 @@ }, { "cell_type": "code", - "execution_count": 224, + "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ - "0.9987894616589735" + "0.9987894616589734" ] }, - "execution_count": 224, + "execution_count": 5, "metadata": {}, "output_type": "execute_result" }, @@ -230,7 +229,7 @@ " return y**2\n", "\n", "def f_i_prime(y):\n", - " return 2*x\n", + " return 2*y\n", "\n", "assert np.allclose(f_i(f(x)), x)\n", "\n", @@ -262,7 +261,7 @@ "\\end{eqnarray}\n", "$$\n", "\n", - "The intuition of taking the modulus is that if the (local) rate of change of $x + dx = y + dy > 0$, than we can ensure ourselves that the rate of amount of probability is increasing, and the total amount of probability is preserved. This is ensured by conditioning $\\frac{dx}{dy} \\geq 0$.\n", + "The intuition of taking the modulus is that if the (local) rate of change of $x + dx = y + dy > 0$, then we can ensure ourselves that the rate of amount of probability is increasing, and the total amount of probability is preserved. This is ensured by conditioning $\\frac{dx}{dy} \\geq 0$.\n", "\n", "![](img/absolute.png)\n", "\n", @@ -272,7 +271,7 @@ "$$ \n", "\\begin{eqnarray}\n", "P(y)&=&P(f^{-1}(y)) \\cdot \\left| \\det \\frac{\\partial{ f(y)^{-1} }}{\\partial{y}} \\right| \\\\\n", - "P(y)&=&P(f^{-1}(y)) \\cdot \\left| \\det \\mathbf{J}_{f^{-1}}(y) \\right| \\\\\n", + "P(y)&=&P(f^{-1}(y)) \\cdot \\left| \\det \\mathbf{J}_{f^{-1}(y)} \\right| \\\\\n", "\\end{eqnarray}\n", "$$\n", "\n", @@ -281,12 +280,29 @@ }, { "cell_type": "code", - "execution_count": 22, + "execution_count": 259, "metadata": {}, "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "(50, 50, 2)\n" + ] + }, { "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAPwAAAD5CAYAAAADZljUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOy9WYxcWXrn9zt3v7EvuXIna2FVdVX1KnWjpda0NGrIhoWxDAgYG7bhJ+t1DBgwYL/Ms18MG/CLBY9hDGSMH2R5xvDYFnrGVm+j7lZtzaquKrKKTDL3JfaIuy/HDycyKpPMJJPMrCqyef9AocjM4L03Is53zrf8v/8npJQUKFDg+YD2ZT9AgQIFvjgUBl+gwHOEwuALFHiOUBh8gQLPEQqDL1DgOUJh8AUKPEcwzupCQggdeAvYkFL+8cNeOzc3J69cuXJWty5QoMABvP322x0p5fxRvzszgwf+EfARUHvUC69cucJbb711hrcuUKDAPoQQ94773Zm49EKIC8C/A/yPZ3G9AgUKfD44qxj+vwX+CyA/o+sVKFDgc8CpDV4I8cfArpTy7Ue87s+EEG8JId7a29s77W0LFCjwBDiLE/53gH8ghLgL/K/AHwgh/uL+F0kp/1xK+S0p5bfm54/MJxQoUOBzxqkNXkr5X0opL0gprwD/PvD/Sin/o1M/WYECBc4cRR2+QIHnCGdZlkNK+TfA35zlNQs82xj4MSsdj1GQUHNNrs6VaZSsL/uxnlsUJ3yBzw0DP+bd1T5xmtMsWcRpzrurfQZ+/GU/2nOLwuALfG5Y6XiULIOSZSCEmP15peN92Y/23KIw+AKfG0ZBgmvqh37mmjqjIPmSnqhAYfAFPjfUXJMgyQ79LEgyaq75JT1RgcLgC3xuuDpXxo9T/DhFSjn789W58pf9aM8tCoMv8LmhUbL4+qUmlqHR92MsQ+Prl5pFlv5LxJmW5QoUuB/K6AsDf1pQnPAFCjxHKAy+QIHnCIXBFyjwHKEw+AIFniMUBl+gwHOEwuALFHiOUJTlCpw5ig65pxfFCV/gTFF0yD3dKAy+wJmi6JB7ulEYfIEzRdEh93SjMPgCZ4qiQ+7pRmHwBc4URYfc043C4AucKYoOuacbpy7LCSEc4MeAPb3eX0op//Fpr1vg2UXRIff04izq8BHwB1LKiRDCBH4qhPi/pZQ/P4NrFyhQ4AxxaoOXUkpgMv2rOf1Pnva6BQoUOHuc1fRYXQjxHrAL/FBK+YuzuG6BAgXOFmdi8FLKTEr5NeAC8NtCiNfvf00xTLJAgS8fZ5qll1IOUJNn/q0jflcMkyxQ4EvGWWTp54FESjkQQrjAHwL/9amf7DlH0YBS4PPAWZzwy8D/J4S4AfwdKob/P8/gus8tigaUAp8XziJLfwP4+hk8S4EpDjagALP/r3S8or5d4FQo+uGfQoyChOYB930cJmz0fXZGIUDh3hd4YhTU2qcQBxtQxmHCze0RXpSyWHMK977AqVCc8E+AzzuhdnWuzLurfQA2+j4CkAjON0uFe1/gVChO+MfEF5FQO9iAsjMKKdsG15eqVB3VYlr0lxd4UhQn/GPii0qoHWxAidN8dh8o+ssLPDkKg39M3J9QA3Xi9j+nmPqge++aOkGS4ccp15ean8v9nmYU3ITTozD4x8R+Qu2LOnH33fuVjkffj6m5JteXPp/+8ic1qC/CEPdDqZJl0CxZBEnGu6v9otf+MVHE8I+JL0PRZd/o/971hc9tgT9pbuKLIgkV4phng8LgHxNfhKLLvhH96ObuF1KCG/gxf/3BNrd2Jqz1fCZRemKD+qIMsRDHPBsULv0T4PNUdHmU63rW7vP+/QZBzFzZJs4kN7fHXF+qUrGNR+YmjspppFnOrzcHZ+rif9Gh1G8qihP+KcPDTszPw33ev1+7bBNnOY6p45g6m4PgRAZ1v0rtOEy4sTHE0vUnfsajPJxCHPNsUBj8U4aHua431gdsDAI+2hpza2dMlstTu8/79zvXcAmTjDBJsXRB14tPZFD3G+KdzgSB5Np85Ylc/OM2NaAQxzwDFC79U4bjXFch4L3VAXMVm5qjE6U5N7fHvLxYIbxPB/5J7ld1TK4v1dgcBHS8iIZrncig7q8iRGnOG+cbM5IQPF7Z8uE8h2bBLjwlCoN/ynBc3V3XBK2yhRAghMCZegErXY+vnKufyf0qtsHFVol25WTGvo+DOY2aaxKn+aHfP06s/UXzHJ43FC79U4bjqgBSKuPcd7ullEiZ05tEp4pjT1N1eJxYu1W2TlR5KCbXfL4QSnT2i8W3vvUt+dZbb33h932WsR/XZrlkcxAwiVIMTXCpXeL3Xl74wp/nYDXhoCfy9UuKAXiwktAqW9zZmxz52vs3loddt4jXTwYhxNtSym8d9bvCpX9GsO96lyyDlxerM0N480LjS3mex4m195/7JP0HXySz8HlEYfBfEE5bP39aDGH/ffzskz0Waw7nm6VDXXxHxdqPG5cXk2s+PxQG/wXgrHjgX7YhHHwfizUHL0q5uT3i+lKNqmMeG2sXpJmnB4XBfwF4FjTqTuKBHHwfddfi4+0xUZLS9WK+frGJpnFkF1/R8ff04Cxkqi8C/xRYAnLgz6WU/91pr/ubhCcpNT1JCHCabreTeCD772McJmwMAs41XIZ+zM44ZKXr8YPXFo+83xcdjhRttMfjLE74FPjPpZTvCCGqwNtCiB9KKT88g2v/RuBxXdonCQH2/02eQ8+L+Hh7xNv3+vzgtUUut8uz1xxlCCsdjzxn1jhTsdV97/dA9t/H5iCYUXBLls6ldpmLrRI9L57d6348Khw5KyMt2mgfjlPX4aWUW1LKd6Z/HgMfAedPe93fJJyUB76/WP/qnTU2+gFZLk9MT9032tWeR5pL5isOhoAffrjNwI8fysPfHATc605IspyaYzIMYn726S7/6sOdQzXz/ffR9WIsXRAmKWGSUXNMVrsqkfck3P6z7BEo2mgfjjMl3gghrqA06othkgdwEnLLwUWvoaEJuLk9Yhyq9s9HtYKOgoSeF01PXrXYa65FlisjeJghTMIUTQgc08CPM7ZHEXkOAg4Z3/77aLgmHS/C1DXON0psDIJTqeqexEhP2jJctNE+HGeWtBNCVID/DfjPpJSjI37/Z8CfAVy6dOmsbvvM4FEu7cFFX3EMkkzimILNQcC5BtzZmxBn2bHubs01+Xh7xHzFmf0sSnNaJXO22I/LI1QcAy9Wp/XeOERIgWFolB3jgQRjo2TxR68vzdzm1a6HQCLhiVV1H5XjeBw3/XHCp+cx1j+rcdEmytj/FynlXx31mmKY5MNx8GTa71yTEnbHITfW+4yjlOuLtWNP0KtzZXRNYxTESCmnFNyMVtmm5prHUla9OOXXmyNurA94516Xu10P1xIs113mKzbw4Al5tKpu7USquked1I+i0z6Om/644dPzNs7rLLL0AvgnwEdSyv/m9I/0bONJT42DJ5PqXKtypzOhM465tlDm2lzlUAfa/Sdoo2Txg9cW+T/e22Cl62FqGgt1Gz9J+eacKn/dXxpb6Uy41/GouibtskWWSQbDgFyCJtTGA0efkE+iqnvcSX1tvsKdvcmhZztYtjtJlePg565rgp4XszkMEEheXKg+8CzPQqn088BZuPS/A/zHwPtCiPemP/uvpJT/1xlc+6nHwYUmBAyDhIWq80jX8/6NYZ9vDmox65rgfMOlWbK42Cyh9tWHj52qTzeYZskiyXJMXUNM73dUacyPVdxdL9nUXIvOOGSSZNzpeHztYhMpJb9a79ObRHztUnMWxx/E49TYjzOynhc/tGz3KDf9/o1kbxxxa3vE6+cbzFftI7+HUZBgaIKb26NZZWK57pyq1fhZwFkMk/wpzNbVc4X7F9qN9QHjKKVdtmeuJzx4ahx10t3Zm3BtvkLPiw8t+pWON1vs+2OnBBxKkO0v5JWOx3zV4XK7MruXH6eH4u+Dz/GzT/ZYrqtTvGwZlNsVLjZL3Ol4mIbg7Xt9WmWLb15uYejakZvX49TYH3ZSPyzHcdymslSv8O5qn3dW+9iGxrU5JbrR92Pqrknfj3Etpd7T9WJ2RxF/9PoSjZJqM35/Y0Ddtag5JlGa8f7GgNdO0Wr8LKBojz0F7o8t01xSd0w2B8HsNUfFs8fFpD1PtZfWXJVoW+l4tMrWLA49auzUwVj2cTPU81VnVgXYxzhMuNRS/PhvX2vz1YtNaq710Lj5pKq6T9r6elSVYz8MUFUN0BDc3B4zDhMmkRL02B1H3Nwek2SSubLFIIgPxekSwWdnlZj+/TcbBbX2FLj/xKrYBnGaMYk+GwR5VHb9uJNure8z8ONjT/6dUfjQhpXHJfh854U2//yddQCqjsk4TBgGCb//6iJ3O96MVbffjlu2dEq2MWuBPQmOC3kel2J7vwdwsAOv6pgkUz2+zUFAxdYZhwlBnNIsWTimTpgoz2t/05IS3jxfZ2sYMgpTKrbOm+frpPkX3y7+RaIw+FPgfgM713C5sd6n4piMgpj3NwZIBG+erx9yv48zzEmYMlexj41x4eEJslbZ4ocf7pDlOa2SRatso2nMXN/7E4mX22X+5BsX+PntLlvDgPmqw++/qph5PS9mbxyx2vNwTJ3a9D0Nw+TIWB6Oz0sc3MAECVGqKginodge3DTPNVxubo+wDY1JKLnYKrHR99E1cYggdLldnm2Q+8o815dqs2v6cUrJ/s12eguDPwXujy11TXChWVI18Z0xFcc8Mrt+XExacYwjXfL9E/xhCbKBH3Nnb8LVdpmeF9HzE4ZhwrevtmdGZ2iCX28M+cmtPb52qcGbFxpcbpePpMNenSvz9r0+hgDbUBp6EsHVdvnITPZReYkffrjD1Xb50AY2X3VmxKPjcJJKx/2bpqYJPt4eY+hwbeGzjazjRbTLNpfbZaqOSlTuX/N5bOgpDP4UOCph9bsvzdMoWUipiC772XW4P0H1YKLrYIJuHwdP8IclyA66uAs1Rb7x45QPN0dcmSuT5ZJbO2McU2euYrPa88lyeSjmvt/Qao6BgJnLe7mttOrX+j7AIYM8KgOf5Tk9L5o9z8HP4DiclGSzb7CTMOVed4ImBFfmylxuqfdad81DBCHX1Ge5kP3P7GnQF/iiURj8KXFcdvlR8fRR/+7q3IO18vtPnePud7CTTcXcGWVL49bOWNF0d8aYhkbLtZhECV0vBgm6Jvi9lxeONLRRmHK1XeaV5c8MdncUsjUMmavYhwzSjzMuNktsDQLe3xjQ9xPCJGWu4vDK8meZ70cl6fY3jv0NahJlGPpnz3n4c2jy1x9sk+TQLpuca7izU/wz5Z3jjfrL1hf4MlAY/Blg4MfcWB/w6e4YieClhQqX2+WHkkmOwklOnePc3ZprTmNufxpzG2yPAu72PNplCwEkSc7buz0Wag7tso0m4L3VPm9eaBx5Ql9tl1npTmahRpBkrHQnXG1XHsgz7I4iVjoev7jTpWKbzFUsdkY5H24OeXW5xtW58gMhyFHvY78+fmtngmPq6Bps9ALevqs2wjcvNA4Z7FLd4dXl2pGe1Gef6fNl1A9DYfCnxMCP+ckne6z3A+qOCUg+3BwyChLeuNB4oK5+Mp33k7bAjnn7Xo8fvLY0jbl7GEJgGxpRqtpYWyWTX60NKDuqgiClOkWZEySZxDY0/vqDbfbGD1YA5qs2YZrNymE111SU26p96LlcU6fiGPzd3S6WLig7OkmaU7ENvnapyc2dMY2SOfsMgGPd9ppr8uvNIY6pk+WS1Z6PkLBYs1ntPhiGFGo6j4fC4E+JlY7HMEhouNZMK14IwcBPDmXXnwT3n4Jbw4D1vs+9TkDdNVhuuGRZzg8/3OEHry0SJhl9L0YMAuYqNmGas1CxcU0DXRP87NMhmoCSpeNFCZMoAakaX5olg493xry31uflpRovLVTRNcG5hjt7DwM/5tb2mL+906VdtmYutBpkYdCdxPS9mHRnTLNssVh1kFGCl2R89eJnJ/PDRC2vzpX5ya1d5io2e+MIIUEiudAskx6YtLO/KT6vybcnRWHwp8QoSEhTScn9rJxjGzrDID5VS+ZRdNEffrBNo2TRKJkIBKtdn8vtEgMv5ocfqt+1SjZCwKd7E6q2SZoJmiWb+arN++tD/CSjapsIoe6hScHH20PmyjZpLqk6BluDgCBOudAs8bsvzR96nrmKrTaLMOHd1RDTEPS9mDiVbA8C0DQcQ7A9DMmznMWaS80xDyXejss37Nf4v3apyWrXp+tFtMsW8zUXXYBj6g8k/Z7X5NuTojD4U6LmmhiGIErz2QkfpRmmrp3Krbw/pu77Ma2KzfYo5OUD2f+NYYCOIDN1rs1VuLmtMvFpJtGQTKKMxbrN3jik7hrkUvLSUhVD17jT8bANQcOwKFkmQZoihGASpVxoqfLivuEcfB7X0vlkd8y9riLnNEsW76z2cWydIMoZBTlpljEIUiQhf/rNi4dO5qPyDarGnzLwY9680CDLJQjFoBOCWR19bxzRmUT86Obuodi/iNNPht9slsEXgKtzZequySCICeKUIE4YBjGmrjEOkyee8X4/TXYSpVxpl0izDC9KkUgkOSM/xjA0WiVz1mVn6gJTE2QSvn99gbpr0vViao7JK0s1KrbB9jCg5hiULBPbNCg7OnXHIs1Uu2iaSj7dHc+e++DzVB2T8lQfP5OS99b7eFFKq2SzXHfQBUiUkV5fqrHccA9RfK/OqWSgQM7yDRJmNf79U/tSq0RnEpHlUmnxxxkfbAxmFYLnpaX1LFGc8KdEo2TxvZfmD2XpL7crpHmObeizuPJxddUeIJYIuNf1qLs2W8OQSqhTL5lcnisTJTn9IObte30qtpoE2yxZrHQ9dE0gpUQTYJka7YrFYs3Bj5Rrv973SbOcT3bGJGnO5jDk9fM1xmGClIL/4Ue3ma9a+HFGo2RxbU415vxqbcDmMKDumiAFug5dL6LuWLy8VCNMM5JMcrldAh4sSS7XXfwoZRQmVGyDy+3ykTX+f/fr5+l5KjzqTCJeP9+Y1fWfl5bWs0Rh8GeARsni915emNWJ94UVjk5K8UgWGRxORqVZTteLGQYJ15dqCODm7piV3Qk3dySTIKXqmnz9YgNLt/nV+oALTZdvX23xi5UuWQ4vzJfp+eqEjdKMMM2IsowXF6qkMudex2NzEFC2DTQhWO/5VF2LMM2Is5xmyWStG/DWSo9hELM9jikZGst1lyDO0IVGGKdM9ITFusP20Geh7rBcdw4RXvZxruHSm8QEic/H2yM+3h5RtQ00TRyq8d/Zm8w2yh/d3D3Ug/CwVuGzwKMYf8+iYk7h0n8OOK5rbXMQnFhl5WCH2M2dEfMVmx98ZZlGySLKJHGS46UqO36u6SKQvLPWZ3sUUrUN6q5JmkveON/guy/M8Y3Lbb59tc181SaT8IevLXKxWeZ8w+GNc40Z9bRdttgeRVyZr6BrgjTL0YWGaxp0/YjNQYgfZ1RtnVzm7I5C2hULXdOougYa6t+8dr7Bb19p8fHOmF9vDtC1w51orbLFW3e73NweYekCgeC9jQHjMD1WvPNgt51qFR4/VEvvpDp4R+FRijjPqmJOccKfAI+7kz9uc8yjXFI/yqjWzKlIgwuM+GRHULUN8hwmUULfS2Yn9bevtkhz+UBXXtUxefN8g74f83svLzAO06lbrcZRf+NSA0PXudOZUHVMOpOYJMspWRqjICGIMuZrNqau4ZoaoyDFT1LGYcJSw2a+YjMMYs43y5yru6R5zuut8pFhTc+LaVYs0ARpBmVbY6nqoGmwOQi4vmQ+cIIfFAlRrcLHa+mdVq76UYo4z6piTmHwj8CTLJwnbY457r4Vx+Dj7THvrQ24vlTFizPGYYxt6NztekzCDNcU6Jpgtefzy7s9vnm5OauRH0dKOddwZ6FHxdYZBgmrPR/XUCQdcpBSMld1+HR3jKmrewghaJSVkkxZGDRdm3MthyDK+J0XFpiv2seKgdxYH1B1TH72yR59P+bqXIWyrZ7HjzPWuj6TMMWLU9WSa+qzE/xxWoVPa5CPktV6VufYFwb/CJx04dzvBRynXrM3juj7MZMoo2LrNEsWrcqDC3BfZ/7j7REfbIzw44y5ssXWtDc9znIGvpKXNnVIJfhhRsURvL8xpO4a/P1Xl7ixPmA45QoYhqDumrx5ocG7q302BwFbw5Cr7TLLdYe9SUSjZLJUc/h4Z8Q4SinbOit7E5I0o+KYpFmOaWjoCNpli91RTJRlGJo4lFA7KAZyfUkZZJrlvLfa59vX5lisOfT9hE92J7y0qDTnRkFCmuVUbZOtQcAoTLnQdHl5qfbYrcKnNchHMfieVYZfYfCPwEkFFI+SrLrfC2iVE358c5e6a84EJzb6Pn/yjQsP3HdzELAzDNjzYhZrDlGasTsKCFKTlxbLDIKEvuchkZRtVRaM0pQly+Z83SHJJO+vD/CjFKZsNaTAjzJurA9YqDo0XJONfsAPP9rmQrPEtXmVKR+HqaoymBppCsMwxpY6jqXq+0mas9IZEyQZryzX+Ie/dWkmmLGP+8VAQG1irWlIc75ZYm8SsdH3We97CCnQBby8VGO+anNrZ0zdUWHMUSf4oxh2p5WrftT1n1WGX2Hwj8BJFs5JvYCeF/P6NIYeRykVx+Riq3zkiKb94RBZBo6tYRs6GgLb1PjO1TmGfkLDNXn7Xp9hFKOjcaVdplG2MAyNSZjxy5Uerqnz964vzIzmV2t9hkFOu2xza2cyjevr5FIt3P1GmoO5BmA6cSbi9u6EvXHEUt1hse5wru5SPyCDfZwYyErH4517Pb56scE4TKg6Jl+72KRs6Xy8M0aT8PJSlZcWldy1BLaGAbd2VC//uYaqBhwk3RzlRTVKKn4fhwnvrfZxLB1b1whTia7BD15bOvQ5P2xE16M67a7NV/j57S5745D5qsN3Xmg/9Vn6wuAfgUft5GrB9NDQqDjGjF9+lPs4ChIc83BhxDG1Iym4+8MhDA3iJEcIQQ64lkGQZFxqlcmlZK0fMAkTklwlsLYHIUhYbsB63yeMM7p+zBvn6yzXXT7aHrE1CPh/0m1KlsbV+SoXmi5ZLslzDjXS1F2LUZjM5K12xhGXWmXiNGdrFDEKxvQmMV6U8Z0X2g+o7l5oltA0wdv3erQqNi8tVlnr+dzaGXN9qcaLC1VeWa7z5sXDLvo4TPCilHGQUHOUMs0vVrrESca3rrQf6kUd9LauL9X46acdoiTjtXM1LjZVB2P9Pgbh/ogux9SZrziMAkVV/tNvXjy2F2JfcOTKXJlXl2uz5zl47acRZ2LwQoj/CfhjYFdK+fpZXPNpwcO42vuLy9J1NAFJJrm5Peb6kmo8ud99fByl1HMNF8fQp1LKY+quxVLNmo6DSnntXI2/+Pk9Ja3sGny4NaY7jqiWLDpezFzVJkxyhCYYeTE/+6SDH6cMgwRbExiWgRel3N2bMPRjrsyVudfzGIUpaZrzi7s9RkHCm+cbNEomv1of8NZqT+ndVx0sXUPXhDJyIVio2Q+cuL/70jwrHY922SbLJe+t9ZlEGUGU8pNbe/zNrV1Kps4ryzWuzVfIc8l81WGj71Mydc43Syo0yHLCOKddsR9Jujnoba31fF5brgNSVRYsVYG40/H4xqUmV+fKD4zoAqi5FnuT6KEJvuc9S/8/A/89amz0bxwaJesQYeYggaZkGVybr3Bze4RjqtbUO50J5xvukfHco5RS9+PJgwm1F+YrrHQ9epOIl5dqM7e7WbboTSJGUaa602o2fT/GEHBnz8M2NGQOm8OQSaTc6CDOMByDkg6a0InSjFGU0J1EuKbBne4ES9f4eHtMmKR8tDWiYhvYhgYZ5Eg6oxjX1rnQdAFJdxxRsgzudb1Z6DAOE26sD7ixNmCx5hAkGbahY5oaXiTZHYdqoyyBY2isTnn5NTdjZxRSdQwqtkEuVT6gXjKw9MPe0XFe1H4uYRJl1By1xLdHIeMwxTY0ND6bmadrgp4fHzmia199+Khy7HOdpZdS/ng6SPI3EseV5vZVXoQQXF+qsTkIGIcJORxZtnuUUurB+1xslnAMnZXuhOW6y1fO1Wfjkm6sD/gX722w1vOpOgZpliNQ/e1+nFO2BOMwJTE0wjQnR3HbxbQRZbHuEiUZQZIhpcSxdD7YHJKmOedbJcI4I5yy45JckqWSxbpDnEssIdANSZpl9PyYqq2DEIcy8IYmuLExRCDRpsMebu951F2TdtXBtozpCCzVUvvBxphm2WRrEOHFKSVbp+cnLNWcqZ5extYw5GKzdOjzPCoJdzCXULGVFh/ImYItSKqOOTuRozSbjeiquRZRmhMmGbah8Xd3eyzXXVolkzSTDPx49r0WWfpH4FkeJnl4fvpn5bRJmN43HkrJK1mGdmQc9yil1KPknXIp2RkFVGyDG+sDtgYBPT8hz6Fs6ax2fUByoVkiSJTvIISgWTbpezG5VLX0880Spqahh4I4zWiUbOZ0wQdTOaosl1i64NMdj1TmeEkGUmAIgWlq6BoYqCGTjqnhhQkV26TuWtRd81AG/ub2iIZr4ccJdzoehhDUbIPtYUiS5YRphqFrdL0Ey9CIs5zdccg4SMnynDSX7ExCTK01HXclaE115X+13j9UYvzeS4fnFB7MuSzXndnG41o6UuZEaT5LkLqmTphk/OC1RX744fZs1HaaST7ZG/GVpQbzFZsozVnt+VxqlWYue5GlfwSklH8O/DnAt771rWdK/HtzELAzCnFNg5pjEKU593oeFVvF06BqzPtu91FjmQ5mjlsVm6vtMoauHZJ8ene1R5jkbPYDTEMjyyWdiRrLfH2xxkrX4+bWiCtzFWquSXeiFFn7fszuJKZdsliqu7iWhqFr7A4jbEtTrbqOSZJm1B2DoZ9g6zqbwwBNCC61S3hxyt2uT5Ip/XimC1/XNdIsVyevoyMRNEsWuiYoWzqTMKFdtnn7Xp+vXazPBkHUHIPtYYqpaby4WGWt67E2CJAoF73nJWSZpFq28aPPDN9PMtplC1MvcXNnRJzlSKmWy9YwpOGq03Y4TNjUP+MVHNapUzmXMMn4yjm1uX66OyaXHBp6GSQZQqjqSc0x2egF1MsWkNMsqbl8fpJRnnEAIgxdPLTeTQcAACAASURBVHCfZ6kPv8jSnwCTMEWDWb+7Y+pESYqUynW/sT7gvdXBsWOZDrrq37zcYqXjKQO51JhlgfeTf6s9n84kRtcFlqahCUGS5dzem5BLxXQbBTHNksm2pav20iQnyXMMXRBnkppj0ihbVG0DL84IkhRdCL5xbY7Vnsd6PyCTOULAtbkyF1plPt0dM1dR1NVcCqq2pliBGli6QZjkXGmVaZRN+n5C1TZplk1MQ6ddsaiXDKSE99b69L2Yuzl0RiHnmw7lqZLuV87V2RkFU88FHFOQ5fl0PFRCzVGeUtk2sAx10vcmETXX4m53gi40xmGGa2lcmasgZX6k7NVR/fH7ZKP97sFgOhpbogRLBHC5XVJtvWnO+bqjyoDjkHK7gm1o7E0iXlz8bDDls9iHXxj8CVBxDLxIqbDux5S5lFQcg0bJmo1lOhjPwfG8669etA65/vubwbX5Cj+5tYtrG5i6ztYwYK5s0q64fLAxRErJSsfHMgXfuTqHY+mEcYZr6yw5DqahMY4SOuOQkmUgUW2xL8xVVQIPyZW5Cv/p773A3Y7HT27tYZs6uVQeSHcSIxBoQpLkkoSMxYrD+brLOExxbEWY+bffaMwShwdLae+u9dnoB9i6IJeQZBnDUI3IutuZUHUt2mWbIM1IM/X5GZpgFCgv6UKrRNU2mavafLozYhTEhNPrZ3nOctPl9t6EpZrNcr2ElBqjMHlA9gqOJtPcfyLXXBPb0ClZBl6cT2P4jJ6nvDMvDhgEKZekZBTE6BoPjJ1+1nBWZbl/BnwfmBNCrAP/WEr5T87i2k8D9ktkfT+e9W8vtCu0Kur0fme1j4ZqTjmqDn9SXrYQgvm64qQnWYYGzFUdwjjl070xr59r0CwZDIKUX631udIu825vgKkpQo5rqZNqfaDkoV5aqLA7Crm5PeIPXl3kqxcbswW7Pc3c3+v6pHlO34vxooQsz7EMjaprkiQ5QZgTliWvnqvzn3z3Cpfb5Zkx/eyTvUN89oqtQp5JlPLSYpVmyeKDjQG3dyaULKULHyY5JdsgkxIN+P3ri6R5zv/+7gYruxNeWqoy8NScvo2Bj2uZaBpcmatgahqaQElsoxJuFduYfZaHKxwBV9uVB6bHHqyr/+jm7qy3YT/BZxs6rqVahFtlCz/K2JuE6JrGD15beupd9kfhrLL0/8FZXOdpxdU5tcgvtkqHEjStsjqdValHud43t0dcX6odqsPvSzopDr0aTXyQQ38w4/vyQpXbuxNsQ2OxahOmkg+3RliGylxL1Fjonh8T7ox5bbmKlGDqBnVX59NwMlXPhSCRLNZdLk5HLF2dK3NjfcDf3u7iRao5ZXsUEEQZkyQly6FWstCFwNZ1dAlLTYfvX5/nckuRVkZBMuuxHwYx28OA99b6XF+q4UWKEJTmOd1JzE8/2WMcJnQnMbWSSdMxaVZVt10UJWwMfP7l+5uACi2GQULPj1npTHhlqcZC1WGx5pDlsFC12R2H6Kh8ycHxUer7yPjLt9fIclUSLJkqPHItfRaz3+8FHPzc1biqMVGSMl9RY7pWuhMut8uca7jPRK/7SVC49CfAcTTKnqcGP+5ryeVSsjcO+Xh7hGsZs4SRoals+HEc+oMz4RxDw7F0/Dhjue5QyiR+kjNXtrE0SDSVhPuty00GQYqhKyLMcl2VrG7telRdA8fUeXVZ3T/Pc+50PN5d7fPJzgQvSnEMjWbZYmsQEus5Ipaca7rYhs4oiOh4EU3XJEklr0wTXbujkL96e53L7RK2pTGJMnZHEefqNluDkEmsNOmkhLfudgnTDNcyyEVElGR4hkYlzchyVePOpaBVdijbOvNVm69davL++hCkZHsY0ChZdP2YubLS2VuquwzDBD2TM9krXRPsjUNWOj5VW6dVUQ1GYZKzPNUguL50NPPxYKa9YhsqC9+dULINWhWLb165+Bth5AdRGPwJcByN8mAd/nzD5RcrHfwoZeAlXG6XGQcJvUnMStfjylyFJMsZR9khDn3dNbmxPsALldvuxxl11+S7L7SZrzpsD0O+faWFbeiMowRDVy7o9ijmUstF1wSdSUSzZE0ltTS8KKFR+kw7fhwqzyDP4YONAZmEqm1ScXTKjsGFVon3N4aEieKqgyTPJZahYRmCSZiyOQj4eGvExztj+kHE9jDCtXRVtsokMkqYrzr8en3IMEwYhikV28TQBRXLJIjTaektwdB9dAFzVZsky9kbZ1ydU6FNydapuyaZhBfmK/x6Y0hnokQ35qs2v/fSPJfbZT7cHPHuap/5qkPVNbAMQc1VYVGjZOJHGaMgQddUyfOoGvn9mfZ9IwflDfxqbfDMKNmcFIXBnwDH0Sh3R9HMJRyFCS8uVFnrBTTLKqkUJil9PybLc+I0OzRySUrJWt/nvdU+760NSPKc5ZqDYxoMvIitYcgfvrYEDGiUTH50c1fNpXNNwkQyDGLeOL/MMIhhKDB1lcB6ebHKR5tDTE2d7PsjoBeqFvd6HlkuyPKMT3ZHeLESj/Qj5R6TS6RUra0CSHLwopR/+vO7LFRt1noe26OQMElVUjCH1a7H+WaZC80Sn+yO+bQznra5QpbnnGuUsAyNIIEwykiyHBA4ho5tZVQTla+4sa6MF1RffJbDWi/ANnUyD0ahEga93C7TmUSHNt9f3OlSdw2iNCOTSuH3460RaSb5xuUGF5vuQ6fo3p/sO41wxtOOwuBPgOOSbhXnszr8OExmevSvTF1p29AZhQmtkkXPTw7psMdpRpRmrPd9RoEaoTwJU15YqNIoqyEMKx2PmmvimDrfv67IIR1Pla++c63NckOd8MMwPZRfaJZUFeDgCOif3+4SxCmCnLdXB+gIDF2VGNe8kIZj0QtUqdGdJq78MKFsO3THIQ3XYG8S03AtoiRHSqXeQ6qx3vdI8ozbex4GAsvQSfOMcZiyOwooO6ZS8tU0JGBpAl2DNM1Z2Z3MhDgNXXCv4zEIUhZrNje3BmwMQwxN8FtX2lyZK/OLle4Do65aZSWy6cUxa12PcZTScC2GUcLuOOLDrSHfv77wwOjqowz5WeXInxSFwZ8Ax9Eo95M5Kx2PHMUzv75UxZi6kftZ5GbJYnsU8qt1lVHf6Pvc2Bhi6hp+nGHrGq2KMqTbu2Nemq9i6aqL7qsXVf24XjL5k6+f5/2NAX6c45oa//LGBmt9n6arJKOuzpV4caHKH72+xDBIZjmHn9/uEqYpe17ER9sjhMzJEMSJJJXg6Aa6oVGxdPwkJUNtYA3XIooVBXcSZRhTLfxmyWIYJIzCFEMTaLp6T2kmmas5GEFCnAXESUZ3LDE0AUKj6uiK1INQmvOawDAEhlC8984kwYszSqbG9jDEi1R4o6EaiL73UsIkTNVAjKpNxTZm38G//niH7UHInY6HrguWag5///oiFUcxF/en6D7KkDcHAX6U4sX5TAG4YhtPPUf+pCgM/gS4OlfmJ5/sHVKOMaZjmPbjvD98dVGRVnK41/OIkpRcShbalak76fDB+pDbexMyCTXbIJUQJTFRkmObOmVbZximbA1DvnapMRsEcZA51q7Y3LvX59PdMVsDJROdpjkIuLOb89tX2wyDhH/2i3skuUST0A8SepOIJJMESY5lGoRpRp5DnGUIKRARpDnoQqNiGwyChEEQY5s6DVctE10XeHGquAeuycsLFXaGIbvjiK4fs1C1qbmqo8/QBNvDgCjNWaq52GaCpYNjqYRdmuVEiRrY0fNi/uDVJc41XH650uWduz2uzldwTB1D00iyHEODv7vbJctzogyuzVeI0mzatKTT92JFMNIkAsHAjxVrb6qzN/DjWRJzH0cJmWwNQwzBjFd/c3vMpVbpSFWiZxGFwR+Dg8QNIZSQ5L5yjB9n9Ceqw2q/zntQby2cuucVR2V7W2WL99eHlCyDxbozFVuISTLleqa5ZBgmZHmGbepYpmC54cxq5vtx5sCP+cu317i+WOW9tQGmoZFKqFuKCSdlwr/+aBtQM91bJRvTECSpRBMxH28NqDkGAz+mZOoIS4DM2RlHVDJJs2zjRyl7k4hcqpbSXEocSxlw3TXZ6PlMgpg3LtSVHJVr8NVLTf7VR9tTEoyJa+m0ShaGBiVbxzI0TE3QnUT4cUiOSghKKXFMjVeXa5yfNsa0yhaubcwaiiSSRtnC0QVrfZ83zzfoevGsZh4lGT/7pKPaayXkUmDqgjDNeOtuj4qtYxs681XnREImV9tlVnseUZpNr5+y0p3MknmPWitPe5KvMPgjcH/i5sbGgEmY8OYFJQx5c3uEpamJqgs1h5JlMAlTfn67y1LdeaBu++5qn1bZ4lN/QpZDydRplk1u70Y0yyZtyyKZUkybZYvrixW+99L8A4tmpeOR5dCqWOyMlRx1JpUyjKHrNEsGP77VIYgz4jTFMAzKlk5lyroLkoySEORSyUK7ls44yLF0DWtKQDFNjaouppNzYmqOwXzFZhAmjAPFNNyZxJQ6Pi/MV/jBK4tUbIM7HY+PNocM/Zgw1Rl6MVGSUbHLVF2TNM/5ZG+CkNAoGViaQZCpJppBkPDR5hBjKpK5VLPJgFcWaqx3fbIsx7ItBCmOpfPthTajUIUUZduYdhtL5qu24tmHCZauMfQT7nY8vn2tzRsXGo8c3z0KEuarNq6lynmjMKFsG5Rs41gDftaSfIXBH4H7EzdpBrUDgoy742jKfvPw4xQvzljrTjAN1SMepzk/+WSPumsipRrs2HRNuuOIzjjCdfQpF93CmjbJXF+s8fVLTTSNGRvsx7d2ubE+pOeFtMo2IKhOGWGmrlpfwzgnynPmbIO9UYQXZ4RJhqUJTJkxCRP2ECzXbVpli91xRNU2SKXEj1LiTPL6uSrr/Zh2xaZRMmfNMpahqZbbKGPkJ5i6YK6qWIeNss1CzZ5pzv3gtUVyKbm1PaY/iUFAo2yx3HAI45wklyzXHfJcognFDKzrBhLIMjkVnFab0OvnG/z41h43t0bUXB3LUE06b5yvc7lVZqHmsIwLwO4oxDI13t8csVizWag7OKbGva43bcbJqbkmddd8ZLPLfq5mv/MRmFGgT7pWnvYkX2HwR+D+rLwmJGu9kHGU4Mcpm30fQxPoQvLz2112JxFzZYtzDZNbO0r8Yr0fMPST6WDEnL+90+Fyu0yS5exOIrww5dp8mYWqg6lrXJ1XVN19N/6nn+zxyc6YQZBg6Tq7o4gozWmUTOqOxcWGyztrfRAatiaIkpwsl5QtHdcS+HGuxC400HSN1Z7PCwsV/DjHixJcy8C1da60bc43XQwjIEozBkGCayi32I9ScqAXhAiUxFbXSyiZyoBBslRzub6k5sb/vZfnCZOMpmswCRX7LUyU239zRwlpJJmk5prMl20Qkr1JxPl2iWvzFWxDZ3sUsDHw+fe+cYHNYcDeOEIX8O2XWriW0gcANbt+pePxb253qDsGq10PDUnFMXFMpQF4db5C3TFZ7fkMg4TvvTT/0PHdR+VqjmrBfdhagadbCKMw+CNwMCu/NQi4tTNmvRfQLJms7E0I0ow4yQmmWXhDqJh5Ps/JpeT9jQHLdZc4y5lEKXujiL1xTJrDq8s1jN0J632P3iTmey8v8N0X5g6dNGqCSUKSSaq28gLsVMeLU7ZHIbd3J/R9NUXVMgwut0og4Np0nJSlC5Ic7u5N8JOcuqERZDnDIFWklDjDsQxsDZBwt+Nxea5Ma1pN2OgHDP0Iy1SsvjSFvh+RpDn1ksV8xWLox9zaHrE3Crnb8TAMjaEXc7Hp8vJijXtdnzudCYMpx93QdSSQZDmmLliqO3yyO6bhmHz7aptRkDAKE/woY6nu8saFBm9cUIKXv1ofECYZLy9Wp6IgHl0v4sOtEZfbJZZqLiXL5Jd3u0RZDlKFPSXLmJX8frndY7Xr8d0X5x8aY0+jg5nKrzjyVUevlX08zUIYhcEfgX3K5SRM+flKB0MILrRdTF1jZc/jSrtEZkLuRwgEFcdE1wR1R80+7/sJy3UHTcDN7TFJBtcXK4reujbgfN3h3JU2wymDbRgkhxbgKEhIspw0k5RtteRMQ+BPUrb6PmGaU7UNWiWbNMupOQYXWyWW6i4fbA3JM8mVRoksz/HjDCElYZKhCXAtnaWaTZxKNkcRzSznP/z2ZYIk46MN1WvvGhqf7EiEkAhdsNmPmYQJUkIexNhjjcm0Eea9tSFBmlN3TAxdI8szxkFCJmFj4KMBhq6xWLXYGoZcapW43CoxDmNA8jsvzbNcd6cTdeDf3N6jdeCz+HR3zMhXVGTHVGWyN87XudvxWK67zFfUoIsXF6q4ptoMbu2Mma86nG+4CKGqBY6h0feSmbTVUTH2SsdjvupwuV2Z/cyP04e658+aEEZh8EdgvxT21x9s40UZSzWb+ZpL2VIcdZmruesvzlfIcghTS/VWy5xhkFG29VnbpmPqNEoGfpTSKtvYus4giCnbJks1G0MIfvjhDn/6zQuzBVhzTZIspzNl3FUcg7JlsDMKSVLFhtM0lXneHYWs7E2Yr9rsjEKiOGOt77PZ96m4BqYm0DU1G67nx0yClCjL0XRBq2JSsZSgxJXpIr+xNiBIMnIN0kTi+wlelBClKmtu5DlrPdWD7pgapqEx8hMqlo6hqdPt15sjKo6BJmASZNhGzqVWmZpr4UUpXpTy4mKVf/jbl+lMIvw4nRmLrmnTfMVUF29jSDwVqri9N2F3HPLVC42Zsm6U5jim8n4GQYxjalxplZirKCFK3YuwDR0pIcoyslyyMQgOCVmeRqfuWRPCKAz+GDRKFkt1h29dbpLmn4lfnK+73Fgf0Cwr5dm1vo9E8JVzDWVwacZvX2uR55KNQUjZMjA0jVu7E4ZBSp5lmKZGo2Ryaa5MydTZm4SzEUyjIMGLU7YGIUp4JmcUJuyOAgZ+zCRKaFcc/Cih48XEaTbNEXSV8KNjTE/QhDCWZDLnXMVlpePhxxmaEFN+fIYmJbau48UZv1jp8sJcmYZrcaFl8NNPOnQnMaauiDGGDkmqJLPRQBNCDXZMMqRUrELTEDRcE4lSx605JqGuEoAfbA5pugZzFZtLc2WyqeaRrgl+vTlAInhpoTKbePvp3oS9ccD2ICCRkrmSzd5EdRwaumC+6tAq26z21HjpvVGAF6X0vIQ4y3h3bcB8xQYJ1xYqTKKU5YbDze3xA0KWp9WpO40Qxhdd0isM/iGoTeWU9heVbSi5p4utEi8sVPh0Z8K5poula0RpznLD5QevLc56xv/6g23WBgGdcUirrL7EzX5ChpyKKzKt92ozAchmyWKt71MrmdRLJh9vjel4yviDqXF5UUKSQ55LcqEIM4auYemaamvVNa7MVdEF3O15ZLnkYsvlbkfNbctygSIDCiqOwWrXp2QbCKE66GxDm+YQMrJcI85VLV7IlFyCkBBnEiFUqJHmks4kxLEM2iXFOwjTDMfS2BlnJLmkLKA7idkZh1ydr2Dpkr/4+T2++8IcVdtkYxDwL97bYBKkZFKiC8HGUEliLVRsSrYOEvqTmI/TjH/0h9e5szfhUqtEz4u41/MZeDHNss2Fsstqz6PnxQRxRrtqca7hUrJ0TF3nfiHLJ9GpOwtD/TJKeoXBPwT7ffD7i2pvEqFr8A++dl4Z9QsPfunArEFjueGwOw6xTZ26Y9EsKZJOYzqjbaPvs1C1SfJ8JgAJqgxYtQ22RiG/dbXF1jBASMmNjSG744iBn5Aj0YVGmuboQjAJEsZBrLyN5Squoc/yAP1EZfqFEJiGQOaSimvhhYni8SdjXj9XJ0wyrs1V+Lu7XUVeSdRk2iST1F2V1Z9MCUgCsHSBPqXb5qgNsevFvHmxgaFr/OJOB1vX0IXAS3Ild23p3NgY8sa5OpmU3FgfcK7hsjkIuLE2JIhSLs65VG2LOFH01pKtvKQwUe67NtWh8+PPCE7nGkrDvuHaWIbGi4s6m30fL0oxdMH1xRprPR+pPShkue+yn9Q9PytDfdK5hafxAgqDfwgOLgBDF7y4WD30YZ+k08qPU6IkI0lz1aNdd+lOIsI0U259q8yvN4d8Zfkz6aSKrXN7NyDLVX/9wE/UYrcNFtAYeDEbYx9HV7FpjkTXdfIsI5E526OIJFdCEUmWYRk6cSYpWTpIlcjzQ1VfHoVKaPLj7TFpDlGWMvBUnO1HGWVH5RMmUYpAqfrk5ARRjq4LkKBNmXOLNYc0UxLXr5yr8zOp+PJJnBFmOQ3XoKQLRoFqn9Wk6mvfHYfc63r4cUImoO+l5FJD0yDMcgZ+wssLVZJM0pkExJkkTnMuNkuzz/hKu8R7axFSqpq+jmC+anN9qcrldplWxeJeb8IkltiGzq2d8QNCJEd9p/dj33MbBAntsjVTOILHr70/6dzC03gBhcE/Ao8Tnx21YzfLFvc6PnGW0yhZLNccarbB7iSiWVaL7WuX1Im4D8Up71G2NFY6Hq5pYGhwse5yY2NIs2IQZjZZnjMKU0qmRtUx6Hs5JcMkSnJubY+ouRZVW6fhWuxMQjVMQlgkMlaGISVplrNQsRmFKb9eHzCMUjRNkCNxTV3RfQ2dJFUNL1mW0yxZlAy1gaW5pFlWG9fFlstCTVGH/Til7piMg5iFuosXJaQ57Ex793UBuRDUXYNbuxPGYYpAwzYglxDGifJINMX3X9mbkOYSP0356oXGASlv1cATZxlX2hUGvhp9VbI0lssutqFxruHy9UtNWmWLf/7OOpamPXKY51HYN75BEDNXtokzyXurAyqOTp5DTv5Yp+9Zzi08KQqDP0Pcv2OPw4QklfSCBJnlTOKMzWHA0vTkefOCUq0d+PEDhI+5qoU3PX3zXP1sGCTowN44Zr5iYeiC1Y6vss8yZ7lhEyc5q/0ADbgyZ6LrgkxKaraBzCVJlmLrar57mu6z2wyEgN1xRJqpk9vSdcIsRSIIEzUi2jF0UimJcomuqfi/Yhu4lop3y7bJlXaZ6nTay9444pd3utRsg7myxZ2OEg15fdnhzQsN/va2EtEchwmGoWGZGjKX6LqqLMg0xbEsKpZOvaQ4+lujkDSTvLfanyVOwyRlbxJzpVUik3Iq8SUZhQmLNXsWah0e5pk9dJjnUdg3vnbZJs5yshx2xiGTSDEsc6k/1ul7kpzBWRN7CoM/Q9y/Y28OAoQA19Qoly3GYUKU5XhxxmvLNcbTVs/NQcDdjodtapiaBlIwX7HZGip9+s1BQMM1cUzViuqYGl851wABaaYSgLqmc6nl8uHmGNcQ6LpOf0rekXW1cAxNY6nmomvK1b6966FpYiYJZZsaSzWHrVGIkJJMqjhdakotJ80lupCUTA1Qz1GxTcq2RpTl2LqGH6VThmCMLjS+/8oCt3bGbI0Caq7B6+dr/P/svVmTpOd15/d73n3LPWuv7uoF6AZAEAtFiuCMpJmRxdFMeEKW5RiHb+xP4G/gL+Bbh8MXjrAd9pU94ZjRjMJWjEbWREiiREokARIklkajl+quvXLPd1998byZqF4ANNANEA32uSCb7Kqst7PyWc45//P7rzQsOq7Ob1zo8ubumCDO0SpoWTpBmqEpMI9SEIKWqbHRtllvy577Rtteoqk32nLgRgiFjm0QpDlFWfHh6Yyua/LKdvseZv1CK7/wpwMJInnUxbNYfJJ/N+N0nuAaKmEm6wILluGjnr6PUjN40sKeJ0Wt/WfA/wCowP9SVdV//yRe92kLyaY7kgMujqw8z+KMC30XV9cwNFV+wCJJi1UU6HtyQq1hSl35wpXm53sTGqbOMEgoypL9ScR620RRBOc7Dm3H4ELfXaK1NKHgGfK0C5ICx5KVdlWB03mFqBTQK078hDgpSMsCIQQqkOYVQZLTU3WmkRw80VWFSZRQVCWaUBEKWJpKw9KI8oKNpolAsNoyOZkluIbC3UmIZ8lBIko4DmMuGlLB953a9fVoFrE3CvnhjSFbbZvvXOgyDaVDTctUabZswjSnqDRe2WxyvufStnWSvKxdaODa4YyilO9lkpeMw4Sqkr59r59rc2soybUfnsyZx/lymOlxF8+9Wvsm+5NjsqKkaRlLg4vPsoHAp6eMT1rY89gLXgihAv8T8H1gD/ixEOJPqqp693Ff+1cZn7UyuuDeXex5jIKEUZgyiTKoKkxN4cPTOZ4pWei2Jrg58Pn+i+v3MNHHYcoPbwyJMnmVt3SVIC7Y7thoqkpeFOQFWIbKNM6oqvpqberM4owf3DhFILjQ9xgECaISOLrC2M9QFOgpJm3bYFJVJFHBJEpQkdZRuiKI0pIoy7AMhZ5rMPQ10rykoKSpSRz0atPiYBJyOktJioK2I+2rh36C48uUoO+ZdGrW3QfHc1RF8PaedNyZhVLHTyVNG989mPLSZoumoxPWFXddU2jZBm9c6nEwiciKCksX9fBSkwsrLkM/XfrzeaYmh3IUwQfHPmVVEcRZ3fqUOK1JmHJpxfvUiblPivuhlxf67j1TlPDkZbVPWtjzJE743wQ+rKrqJoAQ4v8C/jPgqV3wn6cyetYXToQCU1Pp2Do/vTPmUs9lp+dyNI2lp/m5NtuWzkrDrDHOMdeOEoJEqvQapk5eFBzOE7baNklRyPlwIfj2Tpe744CqLHlzdyTNKhomL6z3MBSFk3lMUYKpSuOJIJPf23NMOp5FWVXkZYWrq4Spiiqk57ymCEaBnHJLsxyllEMuHVtnWlf0LVUlynK6nsHRNEFVpGNMUZN1VVXh5klAmpVsdqBhqdwdF6x5BmFWMq295BQhWXmjMGEWZRRVxaUVj9N5wkrD4ngW8t2L3Zrzz1IsM4/l8NJmLcVdbVjYusoPbw7RFVAQWLoqvfgsnawAW5e0mnNdh1GQPtbiuX/xne86NSjzIzebL0JW+yQdbp7Egt8C7p7533vAd+//oqfJTPLzVEZljiz44HiOpas0LZ3jaYRbj6KSl1zouZh1cerKWoNbg4C39yaEacHdUSir6FRYuqz6NkypnHPQEJXANqQ5wuE0rA0YAfNvcgAAIABJREFUpQ3TuY7L3jhgrWUxCBLGYUJZyN64qiqsNDUsVYWqnqbruhRFiZ9Kuq2hKYyCTApmMol3ioqSvmtQlNLOOc0LTFOnLEsatk7bKjjXcxnMY/KqomFKKa1AcDpPqAR0HYOyqtifJhiqACFwDJUbJzOSHCxDYdUzCZKSg0nM5RWX1851+Kvr0qcvLwI8U6Vt67x/PFvKld+43FsaWC7Ann3P5O44wtQkNkxXFRxD2nfP4nxZ6HrY4nnYbW7x+77/hvewVuzTIquFJ7PgHzZQ9IBZ5NNkJvl5KqNNW+ed/SllBccziVU+GEesNSwurri4hnRkcQ0pJNnpufzp29exNJW1psnpPGIeS3yUpgg6DYN5lDMOpff6/jhmtWWgCOn5fnW9gWfqGJoigQ+6KoduclkUdA3Z3zIUhYapcTSJGQYJbVdnrWExjVKEohAkOUkuyMuCjZbLrUGGKkATAj8pMTRI84KmrSFU8HSdrmuy1rLoOiazKKerKTQsg71JgKZURHnJzdOAa8VceroL2G5bDOYR+yPJ/9uoW3UfHM3xLI2yKmnbOj+6OeTmqY9rqOz0HCZhwfvHc7baNr/z/Aqaqshb0hkXmcWNTFMhzgpURSzdb5Jcinc+7qr9sNvcX18/JawJu1khfft+cP2USyveA3CTp81f7kks+D3gLP9nGzh4Aq/7K4vPU9y52Hf5D+8cEdQFOF1VKKmIs5wwLfjW+S7wEVBhFKSst2yqSnq69zyLzbZCx9HpuBZX1jx+cnvEuwcxz695/NbzPYqy4v2jGW1HsuMmNQ7aVFXSrFgO0XimxjCUi7/naoyjnKZjYOmyxXZ7GOIYKq6ukBUwDjMMReAnkl/vWRIxNQ5TNjs2rqYBgn94uc9LGy1sQ+Vnd8eczGLSIsdQZGpwue8R5yVBnHF3HKGr0mZajgpL1n1aVDRtjc22xXuHEgfW1oxl+rPVtnhpo4kQMmdP8oKWpdN1DZr2vXTZxUJbXLVVRfCzOxM6ro4qkKPMac6FvoulK/zWQ+baH3abO5zE3B2HvLTRBBRuD0PyQjri9j3zoend04K5+niUx6PHj4HnhRAXhRAG8F8Bf/IEXvdXFhf7LmGaE6Z5vSDlnz/JSLDtGKw2zPoKG52RzVZMguyB15lFGVttm/WWw4sbTb653cbQFAZBimsoRGnB4TSi65mczBOuHfsoimC9ZZMXBddP54yCjBuDObvDOe8czihKWW13TOk6873LXXRd5WLP4Q+/tcVvXuyx0/MQyMGR71zs8/xaA8/Q0FRFtrdqA4usKFEUwTSUXPsoL0nzincOptw89bm84tFxDdK8ZBQlADQtjfNtKYu9strgyqr0l6uAtaaJZciF3nUNBn7K+a7Lb19Zo18jtDq2ga7Jabmua7LZtplFOYKKu6OQeZwB8rY1i7IH3v/fubLKf/MPLvDCepMCqTRUhGBvHHFrEDC973tA3uYW/nKLOJnFKAIsXWPgJzRMnbZtcHccLTeHW4Ng+fWLW0KaS1HSYihn8hWEYDz2CV9VVS6E+G+BP0O25f63qqreeewn+xXG562MCgHHswhVUWhZOmUFVQmaJh54nfsHcxxdpesanM4SKuDdwxlhWsiZbuQV9bTW0b93OF0KYaaR5N0nWYltqjRtg+fXmqhC1LTWnF7fqJ1h4VzXZhREpAW4pk7DkgYO7xxMOZlGrLVMJqF0dzU12SnwDAmvOPEjfvfqGm/ujvnhjSFNW+Nffvu87EoEGTdO55IZb6m4usGdGl3tpymWplBVgs22w5U1j2tHPm1Hh6rCbtlMo5wX1j12R+HSUOJwElKWFUIoaJp4qG/f2VicstdPfDqOTsPU6DimzOOj9IExZHj4bS4tyqUnfJiWuKZKmn+Uud6f3j1NmKsnZSb5p8CfPonX+qrEo+ZmZx1Lf3x7TMMyMFSFICmI85LzXYftjsM/urp6z/c9bDCnKCu+c6nL0TTCj2URTY6zSo7dzfqKWgFNU2dvEhJlkiCz2jTpuSajKGV/HLLZcTicRFRUmLrK3jhA1PryaSwrekKAo2uoisA1dYIk59TPGPsZtilo6zqeqUq4JIIPj3282ngxyQuGQcov9yec67qsNCx2hz6/2J9QAlkun0sRCk1LoyhKVhsmnqmy1Xak884oZBLJgtrQT/iP7weYuiaHg0qk51zTJEhyLrekR/vNU5+tjv1AJfxsLq4AAz8jTIul/VTTNjj14wcW4cP63CsNk7yU0BBbFwRxQVLkXF6RzID707unCXP1TGn3OWISpry9N+HtvSl745BLfQ9TW1BgM7pdh422zXAuN4OykhN0DxZ7PhrMWWtZzKKMnmtyMk3kKVfjq1sWDIOYoZ9y6qc8t+LgWjoV0HF0WrbBjRMfPykQVFwLM45mCaoQvLzZYBRkmKocRnn/YI6tq7iGwtEsYuSnNUCyBlooYBmCuoUNVcU4yBECikKeWtePJaIrzSs6nsH+NMYzNExNJS9TCcVIc8JKyno7GFBpfPOcy794dWuJ8p7FBaoq0xdLV7h2HCKEoGGoIASuqfE7V1a41PeYxdKEoqR8aP58dqBFVcQSPzaYx7g9j3GYMI9z/ub6KcB9Vfd7b3N/8NoWb+9NmEYZpq7WEFGdyyveMi07u+E8TZirZwv+M8buMOBPfrbPnZE0U/AMhf1xSJAV9BwdR1eZxzlpXjGJEklg6TkPRSudvUVI22kVx6hbcYYcoRUKZGXJrUGAZ2qsNKSmvKxKWrZOlpccTiNUVchnmUTEeUXD0tjsOEzCXJJrM2mxbBsaV9canPgxh5OIwTwhTAsMTcE1dQTQ9xSmcYplqCR1sU1VpJR2HufsjSPCJKPvmSgI3t2f4poarqmx1XYYhwlJUVGUJR1bZ6tl8+JGkwt9qUfY6bm8fr7DX31wwt9cHxIkBR+ezhFCoAlBUlRoKvQ8A8dQ2WjbbGAvC573L/b7B1r8JCfPC1JFMI9LZtGYmwOf7Y7Lc6vuA7+Lh93mFm2/WZTx8pb0BMzLCsdUHkjvnibM1bMFfyY+rdI6CWUeOIkyVj2T3VHIPCnoexpljXw2NAVFUWhZGgoVqqbQcgzujkKGQcLJLOH3X15/oB7wkcWRtHY68SPGfoqpK6w2LVYbFj1POsQO5gm6JqkzwyBFEwq2rrA7iiiqCs9UGQUyf3ctjaqUxaS8KFCQVXhNsbmRBkszSUs38OOMnmuhKGDokhbTMGVBy9I1GqbK0SwmKQpMQyErYRCkZGXJNEpJs5K1piX59boCJWiqYBAkvL034ed7khP40maTvKx4++6EYRizP5bvo6OrWJpKVlasNgyKouTDE5/XznXuWURnf09H05i+Zy4HWixd5u0KCjeHc46mkrir1sKioW9QlNWy8PZxaVvbMbjYf3gv/mFf+7Rgrp4t+DoeRV0njSBKVASGLr3S4qzET3NcQ8pQtzo2YVJwNItxTZWXN1t8WI9xZkXJ/jjCMdV7jCakxVGEJgSaKpnqkzCnKGEUZkzinO+c62IaEr9clCWzKAcEQkCU5QRpRZTmqKpCXsipOCEUyX93dJKiYhrn9D1V0mzKkssr8qp86ifEeUlVwck8olkLWfwaU10JwYqnE6QF4zAlz0ssS2MepWiagq1phJmU785rRJeuKmiqQlqUDIYBWdPmpc0mwyDhf/3rm/yTF9ZoWBo3TgJOphINFiXS777l6ChCFubysuTHuyMEFc+tNphG2T2mkO8fzQjSnO22w/5EFkANVUHUcwoCidbuOToNW2cep1w/mfP6uQ7jmqj7sEX9WdWWT0s//tduwX/cL/hRKq2zSDrBTiI59tp2TQ7GIdMw5UJPtqkUBXZ6Ludi6a66PwlrlxgDo3Z02RtH/O2NARste3lKrXgWoyDlxqnPPMqwNBVFV/ntcyvcHgS4tsrLm20OJhGbHQfHTOk3TYZ+DJWsyFuGtFUyVUUKYjyVNCvRFYVuyyDNJZAyynJ6rsnxNGYaZrQsOZxSVhVpVqCrgqNZyTc3W+RlxeR4zu1hSMeSI7ZxVpIXMr9tWjp+mqMI2TU4mMj+u6Eqsq0X55KblxT4Sc4Hx3PGYcq/+9keDUtnGKTERQFCTuelNXcv8woqBJstmwtdl1GQ8ObumH/147uc7zpcrkUwPdfkZB7zi/0JjqkxCjJUVbDRtMiKgu3tNsezlKL2u0/ynLujkBfWmwjBxy7qp6ny/lniK7ngvygRwyft2o9SaV200k58eTX3TI2mpXPqSzPHFzYay3HMxc/68MTHM+XbnBal5KSnOX/9wSm/9+K6/CDfGeEaGt/cavOWnyKEwDOl40rXNamqihunPq9ud9hoyY3hZJ4wDqRUdDCXrbyqlEqzpChAQY6oKgqG7mKpKlttVyr3dJWbA58PjueUBcthmXGYYTkKRQU7LRtDUxhOE851bPbHMdcHAXlZYWqCoqjIy4JRkFACTq2Rn6UZK66Fa2o0HWlKWSoVYSZNLiZBiqWrZEWtGTBUyrIiyQtUITBVQV6WnPopDUtKe68dzSgruHkasD+OmIQJrqkxjzN0TeHN3TFFVfHNTSkKyvKSNy73+H9+vs/dUUSYFkwi6QVoaII0LwnTnCgrGAUpecHSKXZx1X+aKu+fJZ6E8OaJxhcpYji7awsh7hFRLCqtZ+P+SuvFvouiwAtrTS72JRk2K0t+78U1vne5R1XJnzFZarY75GXJ7iBgbxKwgNqcBnKc885ILqD1pk1ZVvzyYELT0djpObimTpwVvHc4Y+gndD2TJC/4wYcD7o5DWpaGpSkURUlco5oXeuWypMZayf5+GGdMIwmQ9OOctCg5mkTMa196gHmSo6uCrZbDetPk91/eYKVh4xgKa00bVRVoipyscwwJyszyiiAp8AyVpq1jGiqOrmHpglfPtcnyCkORDL2mpZEXkrIzj3M6roGuqWjImkLbNlhtWmj1IMqVtSYvbTaJs4I7o5B39qeoimDFM8iLiusnc4Kk4Ec3h6w0TFY8k1mSM/QTLvQ9doeyOzKNMzxTo1fzBfYnMRf60unmw2MfBUHTko44147m5EW5PGg+7fPwNMZX7oT/Iq9Sn7RrL3zYF/9flBWczmOats5fXjtZ3jTOttJePSexSTdPfUxNXX7f4tYA1AYLov65FR8ezxhHGSuegaWrWLrGStMiSHKmUUqcldw6DUiyghfWm+gKTKKM9dqoYatt15CFkJZtMKsn2dya6jr0E7SaXmvrKh3bYFYXFC1DorAOJtJ73dJk4S3OSgwVdE1lnuZ8Z6eDpgjyosQ2NPqexFY1bRWBgp/krLek/r+sQFcVoqzE0nVe2mhwOs9492CGoUnN/91xRF5VzKKMoZ+RlyUrTQtVQJCVqEqFoctZ+7ysOO+6fOt8G8fQOJlPCZKM/WlUq/ZKFKGQZAW7o4Agzllv2lzoubimRpzlpHnB9ZOElzZaKLWuXhOC9bZF29b5o2+d49YgoOsaCAFCiCWG/NYw4Bubraeq8v5Z4iu34B/1KvV5rv2f1C+9v9IqJJ/xoQv5rD/ZIkU4u0H5cc4fv7nHnXFIlEgsdJTmqIoiuWxJxgCBa8SsNCwQEKU5d4YRpqGgawJHl0YKuqaw1XG4utbk+omPnKdTmIQp0zBlHmU0DZWsBFOXk3hbLYfjWcTr5ztoioIf5wSZlKgO/ZSkkGIcrRYIpXmBZxrYukrPM9jsOPxyf8IsLmiYKvMkI87kjECYVqR5QabJq3+alyR5Qc+zeX6tgQBMXV6rbV0aTYZpwcEkJoxlcXCtadK2NWZRJt1wdI2mLfl/QZrXXn4h37nYq+m+MQCCiqJSsHQJ49gbR1i6HAA6nSf1+y85/7oqRUZF1eDHt4bcHYXoqjTkmEYZb94Zk2RSPLTZsuVmUpWM/GT5WXpaKu+fJb5yC/5RRAyfl+T5abv2/X3xLJcfPD8p8EyVTl3cO3vTeBjH7v2jGbvDoBa4aOyOQgZBxYpnyWtzW46V/mxvXE93KVQVXOjbdD2L/XFEXhSc+gnTMCNfa0AlU4G2rTMKUpqWTpoXOKbOPE5Ya1p4psFzKw3iVF6Z15tOTaZNOZjEDGs1nyYgygs5P24oBHHBLMrl9dbUaFkaRt/j2ok0s7y6bnOu69bts4yOa9QWVrIi7pqSIZeXBXlesd50UETF8TThJ7sjVCFYb5mczCsMTWW1YVGVcHHFwzMlpaflmkzClO2Ww2pTziDsj0N6nommCpRKjvvamiBIcwaBnGR7ebMlZ9HTgtvDgK6jY+gqz63K+frrNfXGM3X8JOPdwwn/+w9SNtoObVu2OQ+mEWkh7bte+xidxNclvnI5/KMMrnxSLv5Jsdi1DU16uxua8rGbxMEkYncUSLfTOsfbHUkJ7dm4P9c7mESMg4S+Z2IZKsfzhIal03UM8rIkzgvysmSn7xFnBdMgYxqm6IpEUtm6SpjmjMKMjZZNUZWczmPujCSW+u9ujnhzd8J7R9OlQq1pSWrsv3h1g+fW3HrgQ2Uapwz8mPcOZ8wT2YIahxnDICUvSuKsIEhyFCGku4SQk3PvHs7RNYU/fG2b/+Jb52hYKjs9uxbXSNorlWC1ZXO+55KV0vJ5EmV4toZrquRFxeEswjNUHFMOBeUFnO86mLoqC32Ghq4JbFOO7H73Yo83nutTlFLya2oK4yDlu5d6vLDepO+axHmBramESc5W2+b6yVxuPqKs9Qsx339pjVe2JerqcBIxjyVkQxGCWZDxNx8O+eHNAe8eTgG43HdpmBpbHZtXttuP8en96sdX7oR/lKvU41RQH3XX9uMchY8spixdJclyyWs7E/ffGoa1QcJO32NvLPvCmqKQFQV+WqDXgx+aovDyVocPj2dEWYGuCdZaNvuTSE6DxQXTKKdh6diGxuE05GQeMw1TKmAcppzMYrZ7Dv/pS5tMo4xrRzMsQ2W7a3NzGOInOZYmBTBJKbXteVmQ5hW6gIiKsqqwDW3J1MuLkkmYcjpPAakfuH4y50LXwzV98rwgyQQ9Tyr+XEPFy0ocXUEFXlhvkBYl0yhjpyenAuO8pO0oJHlF1zUxdAXXUBn6sSw2VhWjMCFMMzbbDmstOeKb5pJg+9q5Drah8sMbA6mNVyrKAWy3HfpewTjIEAKeX/NqM0h5OGy0LN6+M0EIidhOazKPrspKfV5U3DoN2OxIT/kv0vHlqxJfuQUPn74ovwztsmdpBHXuamoqSV5QVhWede9bttig3t6b8M7BhONpjK4qpLWQ5VzX4XSekJUVHVuOoepCkc9rajKHr8Cp23CH0whTVVlr6Jz6GQ1TZRyk7I1jbE3FszUmYU7b1jBUhSgp+JsPh7yy1eLFzRY/uzthFGQURYmCpNqGSU7D0mi5JpausVtDKKhgxdXRNQ3HVImzgigtuDUIOJ7F0rE2kTJhUUHPNaRDrJIgBKw0LHQheG7V4o1LPQZ+wnrL4sbpnI3a6DFwcoa+9IMLkpzDWYQqBA1Tp6hKVEUs4ZRRmjMJU753uV/n/jmXVl2KWtrbc002WoIbpz4Xem49FCOdey/25dctXgsk379TG2aOwpRZlWHpgqqSkM5+jfxq2jrf2GzdI4R6GmbbP098JRf8p8WXUUHdbNtYmrw1zOrWzmrPu8el5GwUZcU3NttcXSv5+9sj3jua0TA0KgVcQ+PSistzKw3++sNTCio8Xa0NE1TSvMQ1NXquwWAeM0tyLnQdep7FwE9IixxFwCTMOJrFqAoUhcz9VQFJUfLukagdVmc0LOlLd/3YZxgk2JqGqavSwspSaVg6aVZQVlI/MAhiOqXOestmHGaczGOSrMDSVGxDYcWz5IZjKKx4Ji3b4GQWyo0rz5glGWGWLyXDC/rPJMw4nScydUilki4rS9quTpDmrLVMGqZOmpdstB0MFbJSAiIXqdyiQHprEFBSUlYq/YaFZ2rcGYWkOdiGrIGMgvSe1O9i36XvmRxOY+Zxhp9IB5+kfr/zsiDOCka+WE4vvr034Wd3JHDzYu9B3f3THl+5HP5R4rPk4p83Fj33c12Hb53vcK7roCg8FIJxtqbQtA2+e7HHlTWPtCiZpwVbHZtXt9u0HJ2XN5u0bR0h5PV6vWWx3bXpuLKa/U9fXuefvrTGStMiy0tOZzG3hvO6c1BRw2cY1g4rtqGx1bIJk5yf7o4kelrU7rBCKvtMQ9CydExN4WiWogjZgmvbOqqiQiWI04JZJK/HQZKhLLv6grgogar2ss/ZGwWc+BmnfkyclcyjlNtnABNd1+BgGvGjWwOSvKSqKom/1lUcU476DuYJ0yhH1yV1duDH7PQ82rbOOEyl1DYr+Ddv3uVfv7nHPM74T15cZ6sjufoHYzlCfDiNqJBOOa+db9/zGWg7Bn/w2hYbLYsoK5iEGQ1L52IN7/jwZE4FvFZvKm/dGXNnGNL3TFSh8MGxf4/u/usQT+UJD198BfWztGXuryk0LJ03Lva5uib7+4vroWMq/P7LGwDLGfqPzBDte3Tcf3tjwI9uDGk48joqTOlQmxeyFWYbCoiK9ZbNessiL0t+fmfCetPiKIqpKjmtlhUlI79gtWGz4skPeVGVxGmBbcjFp6sS9EAQY6pSr68iaNg6TUdn7Kf4cU5elNwdSypt2zGwDGlI8dq5LoYq+PN3j/j+S+vcPJVwi0WlfBxlXFn3qCp492DKWtOisuD9wxkns4i+Z1BUUgDz2vkOOz2Xt/cm7I2jpYvMuwdTZnVd4Jf7Uw6m0XIQaBylbLathxbcdnou//X3LqCKXdJcFvZUgEqO+2Z5xU7PXW7aeVnRtOSGDLIIe2Wt8dQr7Bbx1C74LyMedVP59P7+g69xsS//+2F5Ytsx2GjZ/MPn+kzClKNpTFFkdBydJC04SWPitMJraGy0LUmyaTv4cS5tlCLJpmvYOqoQjCNZ4JuGCWVV4UcFjqVimzqmLs0adQ1UIXBMlbzUiHMJg+xYOqfTWI6mqgpZLmfNiyDBMRRe3m4Tpjn9rsupn/CjG0Mu9F0MVeWVrQ5CCH5+d1ybS1S0avVkkMiayCjIOJnGXFzzuLrexNRU/vzdY4IkI8kqxoHsXDRtOfQTJDNe2mihKoJ5vQmVFQzqPvzHRV5W/OMrK9wZBbxzOEOg8OK6R7dhcfNUWmCd60gH2iQvsOouwSzOvxYKu0U8W/BPIO6vKZzOE24NAzZa1vLvHzbD/Uk6AmnAIAm0uqZg6ipJWpKWFZsdhyCWQp5ZlKGpCrdHIW9c6nPjdI6lhfJ6LRT6DYMLKy5RIo0nNpsWe8TLqnyYFCR5Sd/RcS2dpi3n8U9mcT1plpCk0iLa0jVWGnK8dOAnTGI53Rampdxo4oxf7E0kySZOaVg6W22Hlq2xN4kwNVmg2xuHBGmOqSqUQKthstV2OJrFbLYdBvOEt3ZHrDQtXEOjsGTO33V0yqpCEdBxDDZa0mqqqipO/eShaswFHGPoJ8SmhqXrvHFRwiwLSrZr/fzJLCHKiqWNlHxd0FS+Fgq7RTyVOfxXLc7WFO6OQ24NfS72XM51nIfOAtyvIyjKiv1JxL9+c2/5tX4skdbPrzZomjrzJGeWZFgqnO/aeLaOIuRMfMvS6TgGpq7wynaHvmey1bExNcEkyjmaxSiK9J1XNQXP0DiZxhxPY7K8wNQUSiFwTI2ea9GydTbaVk2gLTB0hXNdl5arY2gqqqKQlRXTMGV/FJAXBTdO/CWv3I8zBDAOUq4fz9FUBU0RZHlFdzH6ahr0PBPX0ri84tFzzSWo8u44oKgVeBVy05nHGR8c+4zChF/sTylKWWMI0pwPT+aczCLevO99PgvHeG7VI0hybpzOycqCipJZmLLZtrF1Fc+ShUJVEVxZayw3tfNd52tTsINfoxP+i261nL269z0Tx5DTXAeT6AHwxSLnn8cZH57MuXY0p2npNGxtuUEIAWUlp8hajsZaw2QSSQ1AkpdcWfMA6eB6eVWOiqqKIMkLNFVwfX+Oqav4ccYszNgVsvo9i1LGkbSEdk2VLK/QVYW2qdGxdVq2tLI+GMe0bJ2drstpkJDnBZqicjiNaVk6LUuve/QBQSZvADdOZ5zrunKxqhq2XlGUcoN7frXBrYHPByc+szDDMSQrr60b5GXF7jBAEYLrJxLB1XMMoprQk+Yld0cB57oO373Y5wc3Brx3OF3m3tMoo2nJTewH10/5rZo1cNbtNStKtroOu6OQN++M2e44XFlrLNt/ixrKrUFAnBV8Y6v1tWrHLeKxTnghxL8UQrwjhCiFEN9+Ug/1pOPLxAgvsMfzWAphsqKk75pMomz5M5u2zuk84drRjKNpTNsxKMqSYSDntssSbg9DZknBe0dzdEVF11Surnt8+2KXF9ebJHnJb17o0PMkVvqD4zl3hgFH05hxlC3lutMoJS1KDEXF0OTJHKYZQZJDBU1HZ7Vp0nR0krzgzjCia2s8v9rg0orHetuiYar15iP75nGek5ew1rJouxpDX7q/dByTKC1IsgIhYBzmCFGRZFI//7svrLLZsuSNopQMvbIqCZKcLC/pugbv7E8wNYVvbrXZ6thkRcHIl/l5w9SZxzlXVxsMw4w/+fk+b+9N8OOMMC1pWDp745C39yb3/C422zbjQMI9L/Vdmqa8He2NAv7ivSN+sT+l637Et/tHV1e/Vqf62XjcE/6XwB8B//MTeJYvLL5MmMGigHcwibB0laKEG6OAtJC+caoieGW7zU93x2i13bOuihr2YNXa74yRLyGU8zBlpWFwud9gEiUYqoptaKw1Teb1OGjb1mlaOrMo5fYooCzhN3Y63DgNMFQNIUpMTZAVFXkpx1yLUlashRCcb9tkZYWpKtyaBpiRVL/NYm3pJltWcnpOV2C1YS7lsVVV0TAEm22XYSB77mUFaVnStOWV/HzXQVUE7xzOuLLW5Mpak1u1V7yaSuPMtqOz1rI4mkhw5y/2J2iqFNucqgld1+LKuiTe3DjxudxzuXXeVh2mAAAfTUlEQVQ8l+YUQqFpa0yjlLWmxfUTn9+5snqP26tnafiJRpgVbHaspdtsXlZc7LncPPVp1UXWr3M81oKvquo9YNnC+KrGlwkzWBTwhoEETXx4GkAFXVfnp7tj/vTtQ77/0hpaLciRHe6Knb6Lo6v86NaA05mEQNq6wsks5u444nzX4fm1Bpttm6oqcS2VGycBO12nVgKWVAjatsE8LghTaZWkqXI2PspLiqzC1qSnu6oIGpZG09Y5mcWMIwmjCJNcjsZWMByHuIbGCxseXVdnHKV4lknHNVlvWoRpyck0wjRUkkLaW4VpThTmTENp3Lg/jrmw4WLpGkXJsrWWtCWeu+caVAh+9wWJ8faTjBuDgIs9lyDNeGdvQl5VvHaug2fq7Po+nqlzMIvpNy02a25/XoBpqwz9lJWmcc/vAiQjYLtjk+TSYMOo24+zOGe1aRGm+VNPs3mU+NJy+F+lmeSXiRFeXAtPZgnvHEzxTA1LV3n3YIamCDZa1hLa+A8u9/ndF1a5djRHFYL9ccibt8dARdsxmAYFaS63BD/NJDU2SYnTQk512QYClrbJO70G++OQSZhyexCSVSWeqTEOUuKiZKtpUVSCrEzZ6shTd+jHEmOtCUAugnlSYKoSdiFExeEs4fVzbVRFmko+v9aUFtgnc04VaRC503U5nScEqYYQoCmCFU9q4jVFZo4tS2Na2zivNkxWGyZ+nOHVVss/vDHgwxOftm0sC5qWkaArgjujkK6rEyYFmgqTIOW5FZdhkGPrimTyIQt833uud8/v4qxK7+p6kw+O58tN0qshnV8Hms2jxKcueCHE/wesP+Sv/ruqqv7do/6gX6WZ5Jchxb2/KPjG5R7Hs4i+Z/KzvSmaItA0wUpTwiR2eiZv3R3xz1/e5MqaxzsHM35ye4SqQssyGYcZVQUNWyPJZX6flRV9y+S//LYs/i3qEmc3sq5rIpBpw+7IZzRPKRHYqlzICoLtjhybDdKcKF0s7gULXiWNUmk2CbRsk4al8Z2LPVmVn0ToqmAW56y3LIpKavUVUdV2UgrbHXtJifWTgqEfstm2aVgaozBlEqW8stUizkr2xyEd1+T9oxmjMIUSVtsmYSYLlOstizgtluCOJC8wKzn+ut6yKauAeZSDqEiynHNd9x4BzqKYuvgMqIrANVRmUUqFYKfXAL4eNJtHiU9d8FVV/d6X8SBfZHzRMIOH9dVvnvo8v9Zk6CcM5jGrDYuOZ6IiUHXBWkOy6W4PAk7nMeMwpeOZ+GnOMEjlFJuuEmUlSlWx3XXZalmsNq3lcz9sI1MUeHGziQBsU+rmq1JaWs2ijEoISuRU3DySxTtDVUgzAYqc7HNqnX/PNdhqW1zoyVn4rmviWTrnus7y53VdA0URHE5jBBXf2Gyy03P5RV1Ma9s6ZSHJsyrwja02Oz2HvKzoegZ/+K1tfnRjSFaUdGvrakfXOIoiqgpsQ+VgGiOEoNWVtYo4y7Fq22lFUXEthZWGjWepfP+lBxHg938GHFPeNC723Ht0+1+XXvsnxa9NW+6LlOJ+XFFQVQos3eb51SZlWaAimCYprq7yw5sDoqzge5f6bHds/u3PDkiznL5nMYtymWtSIIRCzzO41HcZhSnPrTXu+zc9uJHdGgSkeYkQMqe3dJVLqwkH04iDccjbe1PSogSgKKFQwNWkLZQfZzi61Nmvt2zWmhaXVjyKouTUT/juxS7vHsw4nUtazxuXe8tx1EW8dWfMSsOi65q8vTfhFwdTRAWbPZuraw0UBV4995Hu/fYg4MWNJn6S89bdMSezhDgtSPKCaVSx0bToNwzCLOd4HvPyVhNDU9FVhWmYoSga39ppLwGiHxeLz8Dr5z9i23+daDaPEo+14IUQ/znwPwIrwP8rhPhZVVW//0Se7CsSj9K//7iiYJwVvH6+Q5yV/Nk7h3hmCVVFXME0ylhv2vzbn+0tZaZtxyDJCjqujp4I/CjHsxV2eg4tW06VzeNsydjrugajIH3g2S72qQuHKX3XIM7ktN3Lmy3evD3EUOUpXtRW1YqoUBSVnqtzexRwPMu40Ld5/VybnifJsW6tzBv4CRf6Li9uNJc3mfur27MoQ1PkGOtbd8fYmsp622IaZfztjQENS7tHl3C2mv76uQ7XT+aczGLmacGLGw0urzRwDY04K7g98DE0lVe3PzqNH+ZI82nxdaTZPEo8bpX+j4E/fkLP8pWLR0VpfZqW/p9/c4OXNpv8qx/fZegndFydDSo5c66q0rLK1Lkz9NlqO7yw3uTdgykKgtfPt+l5MldebZhLxt7pPOGvrp3w8lablYb5wLMtCod3JxFFUeGYsnAYpSWuJX3gJmGKa0ooZZhmNG2Viz0XXRN850IPW9d4bvUjccrth9xk/Djnz355xHrLWm46QsDb+1Pm9cLXdYWbpyGqAg3ToCwrJlG6fN6zqYlnaryw3qTrGrx7MONcx8GqN884K2jZOnl+bwno16Xg9iTi1+ZK/3niUfv3j1IU3Om5fPdil44j3Uz/jx/eRleks+o8yvAsjX5pEqQ5v3Ghx3bHQVcFl1Y8mrbESi2850ASb1r1KOlq03rg2dqOwRuXe/yff7eLn+S8cxBx7XBOkOa0XZ2sqAjTAttQ0BWVspLPfmnFY60pZbkCyZU735PYMc/S7vFSn8cZu0OfrGR54r91Z0yclQgki75hSpMLP8lomjqzJOXmMOXF9QZluXjezgOpyW8/v0LL1nn/cM7NgQ/Ifr6myeLn2fh1Kbg9iXi24D8hHrV//6hFwbM3AQGoqvywdhsGYZqjKQpNR+PFjcYS/rB4jb+8dnLPYvMTib+aJx/x9M4+2yRM+Yv3jtgdhgxqGKZlKCA05lEBtqBhqvhpgatrtB2Ni32X812Xq+uyTrA/DjmexTy31ljWBqJMutoeTCLeP5pRATtdZ8kW9OOcv75+giIEt4chmiKwdQVTE0zjFE0Dx5DX791RQJwXZ24lxjKF+vndCUGaMw4Sttu2/LfGGWM/pdmXN46vEz76y4pnC/4T4rP07x8lJzx7E1hrmAwDCNOClqXjNFT8RLLtDO1Bh9L7n8WrnVcWPeyzz7ZIRW4PQmxDxTFUgjRn1XXqcVtpB2XpKq6l89xKA1MTlJXkyDfq1zzfc3lurbGkzlzsww+un7I3Dmla0rkWAX5cSLAl8NbdMdeO5qx4NYo6zrB1g2mck5cVqqrw8maLrmsyDZN7GIH3p1B3x2FtzimYJ0XdIXCX79GvW8HtScSzBf8J8aT792dvAjt9h4Efc3WtwVrT4ngWMwlTXj3XfaRn6TgG++OQc11XYprPPNsiFTE0hTSQRu+OoRGnBT3XrMUrAlNXeG27zT95YU2y7IYBau388nFGHIoiGIcZe+OYJMvZ6bmYusIPbwyYRCkfHMxkqysryOKajBv7hFnJ+Y7D8ysNLF3l+vGMSZjWltbyVvJv3rzL7jBEV1XOd22CtGCtaWFoUjADUso7rlV8Z+PrzKF7kiGq6kvVwABSePOTn/zkS/+5nyce54P0ad+7Owz40Y0hd0YBYVLw+k6HCz13uXjvLw7e/3ofV6X/y2sndByDt+6O+fubQ/wko6qkbfNi0q2iYrVt889eWq9n2nMurXjL1xNCdhJWG9Zys7s18Pm7WyP5GraOrircGviSGy8Eu4NAQjh1lbWmzTzOOJknKKLi9XNdrm40OJhExHnJqmfSdWXLsOVIBd1Pd0fS8VVIt5iiKLnYd1FVhd/YkRvhoiJ/dsGfvRmc3ZgfdwDmad1EhBA/rarqocNsz074T4nP2755lAr/Ts9lp+c+oJj7uOLgw57l/h44IKvkexPCtMAxVfwklzm/IZnxpqGy3bH5xpZ0h3XMj1KIlq1zaxDw5p0xpqbQc83lzP47hzPyvMI15FBQmEmX22GQUFFBjeBOC9gbh7RsOYmXpAVR7WFnaQoNU4p34qzg0orHzYHP3WFEzzNRhMDQVOK8ZOgn/O2NAc+tNpmtpsuN6f4b1hcxHPV5zU6+6vFswT9mPI799CIe5l6zKJjBR+DMRzltJmHKNMqYJzktS+fqWhMh5jSsgqatsdl2eGW79VCRytkPuQIkacFfvHdMzzOJ0pw8L+m4Omkhp+8MVeEgSimqijcu9Tmdx/z9rTFpklJVMItrS2nHoOcaeKbGyE+ZxSm66rHTc2lYss3mJykvbrTYHcpC3mAekxcFbdeg5xn8dHfEa3U1/0F60AgFZckGbFj6Y7fqntlFP4sH4vPYT9+tzSnOLtyzBbnFHL0A1poWaV7y7395yGCeoGsqXUdaVk/qPHYBvVxsBkfTmL5n8uq29JJPC5VXttuc7zr8zpXVT/y3/Nkvj5hE0s755iDgcBLjmhqqIkiLikmYsd622Ol7DOYxQZJTVdD1zJq4o7I7CpnV1NmkLOmrUlu/2rAoK7iw4kLFMicH0DSBZxqoimC1afGjmyOmYYpnaby20uR7l/sPFdcs3n9DVVEESwfYq+sN1Pq28Xnj62oX/WzBP0Z80inwsAr/6TzhcCoHas5uEJdWPG6eyl7z/jhEICmyilLw5p0xP9+b0HcMvrXTJclL7oxCzncdbg2Cpapusem8fzSTVlOWRlnJav5GyyIvP75W8xEKKluO9B5NI+msSsWHJz6OoTCLM8pxxXbb5nzXZRaly0Ul/2xwue9x43hOzzVouTq6UBj7Cec7DpoiaNkyX//53TFZUaKrSq03cDmYRkyCBF2BtaZJw5IMu3ntC3D/Ylu8/5dWPOm6owtMTeHmwGerbT9Wq+7LnLD8MuMZ0+4xYkFUORu2rjKLsod65L17OKUsK947nPHB8XzJPB8F6ZKJdzyLEQIQFbqqEqc5FBLSGKZFbTGtMgoSZlEmRz9LuHY04//+6R3+4/sn/NkvD/n7WyNURZAVJb/Yn/BJyIKPUFAyFQmSnFmUk9X+c7Oagnux79KwdfbHkfS2r+APXtviD17bIq8q7owCjmYRl1Y8FFVBUaSbraoq3B4GnO85vLLdxjFVECDkf+CYGm9c7tGyNahz+LWmxctbbTquwcEkeuhiW7z/DUvn6noTXVXkzSIvHzvXfhSPw6cxnp3wjxELVNX+JOTOSF7VVxsmL2w0aTsGl1Y8fnRjyOk8RiiCa0dzep5B25bk13mccWWtQZwV9xTk3jmY0nHkHH2UVbi2Rl5UnM4TXFPik0/9hOfWZNX79iBgbxwu7ZfjLOfuKKDjGlzouVR8/Go/mwPHecG1E58V16RpaUyTnCgtubrRQFMUFAWmcUauqzRtjT/61vZyUX3/pXX+/N1jVEXhfNfGMaX2ve+ZOKaKEHJk99YgYLVhcaHnLZ8hTHNGQcrllQbf3unhJznXjmaoChiqYBCk9DzjgRP77CksF73+uXT1D4tndtHP4oHougb/4ZeHzOKMlq2TFfDBiU/T1tkdBtw89bnQd9nu2Pz7Xx6S5mWNl6o4mkast2xuDQK+sdVavubFvstff3BC3zOpKoWsKJlFGUVZ8e7hlDjLpCmlqtB1DT44mjMKUrKiwjV1grSk37SlhVJtL/3Gpd7ySn8231+03hY58Mk8x1blUI2qCpQawXX9aIaqquz0HL57oYci5BDN2RgFKd/catG2dW6e+my2pIGEpgrWmtIFd/FzPy43PjtEc3W9ycEkYhAktG3joSf2F805+DoO2Dy70j9GjIKUjmfQdAzyUsIjXlxvkpcVP7oxXOb3h9MYVVU413HqhSIr3IN5wn5tZ/yX106WH97XzncoKziex2RFSc8zWWmYdUHPh6ri9XMdbp76KIo0YZxHGeMwqRHXGV3XYLNl0/MsNFW5R4G3gHneGYbsjSNWG9L4Mcxyzvfk4mw7Buc7DpahMA4z2o5OWcHdccgHxz5BUixhkXAvMHIa5VT1Yp+EGXFWcLHncjCJOJrG/PDmkGtHs6U6b3FdP3uN9kyNc11Jll1M1d0fX4bl2Nctnp3wnyMWp+TfXD9lHKZc7Hu4pswvq6piFmVMkpgXN2Ql2q9bZHlZ0SplNTpKK6ZJxnbbXk7AnS3iFWXF/jjiQs/hYBLxzsGMpq2z3ZE45UsrHmGaU5YJLcfg7f0peSF75AKFeZSRF9U9Rgr3FxnzsqJl6czijKvrTUZBysBPKEo5yZfkJYaiLP3tJ2FKkhWsNCxWGgY/uzNetvfuPZ09jqYJkyijYWlcXW8QpQWH05jVhsm14xl3RiHvHkx59VyHpq0tr8uf9Rr9uKfw0yqu+bzx7IT/jHH2lFxrWoDg+olPkH7EjNc0wUpDGhiCrJQ3bQ0/WbDnHDbbFo6m8tJGa8lvWyzGRRFvEmXcPPWxdJUrax5X1hqoQhCk8nVtXUUIOJpGdF2dnqtjmyogWfOTOLvHSOF+hPbeOOTOOODUT2hYOt/casvXLgVDP2Ue5YzChKal8sHxnKyQ1YALPRfXNOh65tJkcXE6n8zipRNMUZS8vCVtoW4NfVYbJkM/Yavt0HcNkqLkrbtygztrs/VloaK/THz5VyWenfCfMc6eklsdh1M/YX8csjcOON9xmcYZ2x2bV7bby1bbRsti4Etjh56rc+onqDWKaqVh3vP6i3y27RgSAmmotByT3aFPXlTM45ybA5/9cYRrqjRtnZ2ey2rD4s44JEhy+g2TnZ7D6+e79/TeF0XGO6MAS1e50HN573BKEBc14lqmCwBZUWCqKpauYekakyhGEYKdnoRfxlnBlTWPWe0YuyhS/vm7xxRlyTc2G6RFybWjOa+db7NRO9zKLoNG1zXZ6Uk/ulGQLhV+X+ZJ+3UV13xSPFvwnzHOFp0als5r5zq4hsr7x3PiRs43NpvLa+7iQzzwEwxNIS8riqrixY3GsmL9Sb1ez9II0pw4K+h5Jm/ujrg5COi7Jk1LZRpl7I1CXt5qs9a0eGFDFv8WJ+xm277n2S/23SUPf0Ft3ar749eOZxxOY17dbvOLgylX1hpMI2lWEacl602TUz/hVaONrkr4o6oIHPOjS+KicHf237Oomjcsnb+/NWTFs5Z/l+QlXUfnYBKxN5Z+eHleoWmCvXHIb9cOMvfHk7qGf13FNZ8Uzxb8Z4z7VXGLYZDXtjsPFJcWuKlJmPKNzfY9lWT49CrzZtvG0tTaL12y1VcbJq6pYWgqL2663D6dczSNsHS58ExNEllVRXYR3rozvmdhbLQswiRfoq0l712KWioEtqZAJVn2TUuwPwkJkoLLq/J7v3W+87EV8U9aQK+ek+YbC4FOksse//muw+7IJ8kriae2pS323jji7b3JA+rAJ6lx/7qKaz4pnuXwnzHO5qrvH07x4wxdVeh75kPzv/uNIxd/vjUI/v/2ziU2rrOK478zb3vGzjhO7CS246ZKZVpoWlpEhYRUoRJUQaGV6AYWVHSBWCDRXSmR2KBKVJUQaxZILCrBgiKkAupDIMSCINoKQtOmKA9C82rsJH5MPBPPjA+LmWuPx/O+9+Y6956fFCnj8UzOncz/ft93vvP9T8cs82K9G+vJS0uUylXumciSScU5uHuYz9w1zux4jmwqwex4lpVSmYO7syRiwnyhREXhkUPjnJ0vbFufjmQSHBzP8vDsGHP7RhnJJDe+5PdM5FgulRnLplharZ10S8aEwxM58pkkB8ez3KpU22bEHQE5OYJ3zl/nxIVa0U9+OMXR+yap1GcfifoN6dy1m7x/cZnl1TJV1ZoXfTLOrkyS01dXNt7bEfqr737ExRtFquu67fMc9P8ybMU1nbARvk8ckb7+3hXK6zCeTW4c2GjVvaTbtLFVlrlxFHt4doyTl5d47cRFbqyukYxnoKFKtrqu3D+dZ3cuRSIuHJ4c2WiK2Gp9equyOcNoNVIvF8vEYsK5hQIxWSeTTrJ3JMXeXGbjBtF8Ft3h0J7sFoOMVDzGUqnMUrHM4uoas+NZnn54hnMLN7m0WOTyUpFD4zmuLBapVKucX7jJ7J4s2VSCWj8e2fZ5xIgRk1ploRPPoNPwsBbXdMKta+3LwFeBNeAM8G1VXez8qjuf/HCKfbsy3Lt/dEubrVZfvEGmjY1iXSmVScbizE1uusR+cGmRT+wfpbquLBXLPPXQ9LZjssvFxY5Ouo1f8n27clvWxKNDSZZWK4jUbhQTI2kO5Ida1rM3fy6jQ0lyxdoR2Vw6xgPT+VqWvsFrr7nL7sHxYc5cLZBOKPPLJeK7MiyXytx3YNe2zyOXSVCuKpmkcGmxyNy+pKtpeBiLazrhdoR/E3hBVSsi8hLwAvC8+7B2Pr0KeZBqsMZZgdOUMp2IgZQ5+sn9vPPfG5y6vMyRmTG+cO/kFrE7Ca0z8yuk4jVTSseyqtFJ1/mSt1oTr65VeOTu3VtMM6GWgOsmLFU4MpXfciN0XGraXePhiREKtyqsFCssFErsHUkzPTa80UGm8XcP5If48MoK6USMlVI5Uk0kvMCtTfUbDQ+PA0+7C+fOoVchDzJtbLyZFG5VGc0kuFWp1k++DfGVIxlurK7x6Fz7hNbc5CgnLi5x4sIN7p/K920e0W3q30vsDq1uhM118A/OjHF2ocCeSnpbb/aa7fUilSrk0nGm8kNcLZRYh5b+f0Z7vFzDPwv8ut2TQTaT9IN+hNzvtLHxZpJNxep90DbdbdpNYZvF+8B0nrMLBU59vMJDB8daxtetiUa/69tuN0JnBtK4ht87kiYeE6byQy0NLpaLZQqlMqOZWjOOMwsFpseG2m7bGe3xpJmkiBwDKsAr7d4nyGaSfuHX+m97H7RKT33QmsU7kklyZCrf0vTRoVsTjX6vr9ONsHEGMjM2TCYR51zd4eZA/fx6s4DPLdzcaFt1abFYs+pKJyLRy90PXDeTFJFngCeAxzQIR8yQMkgftEEShH6cOGt3o2iegUyMZshlEttMKRtxbmIiwty+zfMKYS6O8RO3WfrHqSXpHlXVVW9CMpppFpAzUm7vKde/eP3ammpVDdeLd1+n5h0OYS+O8RNXNtUichpIA9fqPzquqt/t9ro7yaY6aFpZU5+dL7S1ZN4Jp7/a2UbHY7KR+W/07sumExvtrFqbVHpvQR1mfLOpVtXDbl5vdKbVltmb71/h0Hiu7YGPVi2bbrfwe8n8O959ijA1NtzFmjtaxTF+YpV2O5hWwqmuw/Wbt5gY3TyE0lzwE7Snei+Z/4+XS0yOZpgaG96oE2hXMRe14hg/McHvYFoJZ/dwkutNomhe03px7NPN0qCXzP9Kqcz/rtfcc3LpmlOOW2tpozt2eGYH4winkd3ZNPFYrOOBj05uur3g1hii26GUxr31VBzWKuv868IiV1dKoT64shMwwe9gWgknFoOj90129HFrvFE4ybG/nb3GlaVST6LtdMKvF7p5zTl760emx0gl4qxV121v/TZhU/odTKeEVauecg7O9lyhVOH89ZvEgGSMjSO83dbyXhhDdFp32956cJjgdzhuqt1ef+9KzfU2m+54hLcZv/e+bW89OGxKH1KcI7yfu3vPxrlx6G0t77cxRBSNJ3YKNsKHmEFHUi/2vjtl+W1vPThM8CHGTZ28m73vXuoAbG89GEzwISaokdQP++edUDIcBkzwISeIkdRr++egKwfDhCXtDM9pVTDkJgvvti7A2MQEb3iO11l4t5WDxiYmeMNzvO7q6vWMIcrYGt7wBS9zB373gY8SNsIbOx7rA+8dNsIbdwS2b+8NNsIbRoQwwRtGhDDBG0aEcCV4EfmxiJwQkX+KyBsicsCrwAzD8B63I/zLqnpEVR8EXgN+5EFMhtEXTuntXz682pcVVxRxJXhVXW54mGVL53LD8B+3/ntRw/W2nIi8CHwLWAK+4Doiw+gDP07mhZmuI7yIvCUi77X48ySAqh5T1RlqjSS/1+F9viMib4vI2/Pz895dgRFprM6+P1y1mtryRiKzwO9V9VPdftdaTRle4UznG119VtcqHRtUhp1OrabcZunvaXj4NeCUm/czjH4xf7z+cLuG/4mIzAHrwHmgayNJw/AS88frD7fNJL/uVSCGMShWZ987VmlnGBHCBG8YEcIEbxgRwgRvGBHCBG8YEcKzwpu+/lGReWrbeHcye4CFoIO4Ddh13nnMqureVk8EIvgwICJvt6tmChN2neHCpvSGESFM8IYRIUzwg/PzoAO4Tdh1hghbwxtGhLAR3jAihAl+QETkZRE5VTfx/K2I5IOOyUtE5HER+VBETovID4KOxy9EZEZE/iwiH4jISRH5ftAx+YlN6QdERL4E/ElVKyLyEoCqPh9wWJ4gInHgP8BR4ALwD+Abqvp+oIH5gIjsB/ar6rsiMgK8AzwVxmsFG+EHRlXfUNVK/eFxYDrIeDzms8BpVT2rqmvAr4AnA47JF1T1sqq+W//7CvABMBVsVP5hgveGZ4E/Bh2Eh0wBHzU8vkCIReAgIncBnwb+Hmwk/mHNJDsgIm8B+1o8dUxVf1f/nWNAhZqJZ1iQFj8L9dpPRHLAb4DnmuzXQ4UJvgOq+sVOz4vIM8ATwGMarmTIBWCm4fE0cCmgWHxHRJLUxP6Kqr4adDx+Ykm7ARGRx4GfAo+qaqh8t0UkQS1p9xhwkVrS7puqejLQwHxARAT4JXBdVZ8LOh6/McEPiIicBtLAtfqPjqtqaEw8ReTLwM+AOPALVX0x4JB8QUQ+D/wV+Dc1M1aAH6rqH4KLyj9M8IYRISxLbxgRwgRvGBHCBG8YEcIEbxgRwgRvGBHCBG8YEcIEbxgRwgRvGBHi/2Ar9L3OMGT6AAAAAElFTkSuQmCC\n", + "text/plain": [ + "" + ] + }, + "execution_count": 259, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD8CAYAAAB0IB+mAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAbL0lEQVR4nO3d/48cd53n8efLyTiOsOOc5HhtxfkicdFFBxKQs0y+SCsfghUxFvklkbISoM0vViJOAt0ipPADK+4PQMDuJb4RrJYIWIQOkosiJ1pLtwj4IUGx1yFkbZ3CriFzts+QFXasRLZC3vfDdJN2u2e6uupTVZ+qej2kkWema6Y/U55+1qc/XdOtiMDMzPpvQ9sDMDOzZjj4ZmYD4eCbmQ2Eg29mNhAOvpnZQDj4ZmYDUSj4kk5KelnSMUkvzrhckr4h6VVJv5B0R/qhmplZFVcvsO1/jojfrXHZvcBto7cPA4+P/jUzs0ykWtK5D3giVj0PXC9pZ6LvbWZmCRSd4QfwD5IC+B8RsTx1+Y3AaxMfr4w+d3pyI0kHgAMAV3H1f3rP1f+u1KDNzIbq/Nu//V1E3FDma4sG/56IOCVpO3BY0omI+MnE5ZrxNVc8Z8PoQLEMsHVpe9y97YGFB2xmNmTPnXns12W/ttCSTkScGv17FngS2DO1yQpw08THu4BTZQdlZmbpzQ2+pPdI2jJ+H/gz4JdTmz0NfGZ0ts6dwLmIOI2ZmWWjyJLOnwBPShpv/72IeE7SwwARcRA4BOwDXgXeBB6qZ7hmZlbW3OBHxL8AH5jx+YMT7wfw2bRDMzOzlPyXtmZmA+Hgm5kNhINvZjYQDr6Z2UA4+GZmA+Hgm5kNhINvZjYQDr6Z2UA4+GZmA+Hgm5kNhINvZjYQDr6Z2UA4+GZmA+Hgm5kNhINvZjYQDr6Z2UA4+GZmA+Hgm5kNROHgS7pK0j9JembGZXslnZN0bPT25bTDNDOzqoq8iPnY54DjwHVrXP7TiNhffUhmZlaHQjN8SbuATwDfrHc4ZmZWl6JLOl8Dvgi8s842d0l6SdKzkt5XfWhmZpbS3OBL2g+cjYgj62x2FLglIj4A/DXw1Brf64CkFyW9eOmdt0oN2MzMyikyw78H+KSkk8D3gY9I+s7kBhFxPiIujN4/BCxJ2jb9jSJiOSJ2R8TujRuurT56MzMrbO6DthHxKPAorJ6NA3whIj41uY2kHcD/i4iQtIfVA8nr6YdrVs2l23cl+14bT6wk+15mTVjkLJ3LSHoYICIOAvcDj0h6G3gLeDAiIs0QzdaWMuB1X7cPENY2tdXlrUvb4+5tD7Ry3dYtbUa9ST4gWBHPnXnsSETsLvO1pWf4ZqkNJexrmfXz+yBgKTn41oqhx70oHwQsJQffGuHApzO9L30AsKIcfKuFA98cHwCsKAffknHk8+ADgK3FwbfSHPhumPx/cvyHzcG3hTn03eX4D5uDb4U48v0z/j91+IfDwbc1OfLD4Fn/cDj4dhlHftgc/35z8A1w6O1KXvLpHwd/4Bx6m8fh7w8Hf6D6Hvpz772m1evf+quLrV5/HRz+7nPwB6YvoW876POsN76uHwwc/u5y8Aeiq6HPPexlrPUzde1A4PB3j4Pfc10LfR8DX9T0z96VA4DD3x0Ofk91JfRDDvw8XTsAOPz5c/B7JvfQO/DldeUA4PDny8HvkVxj78jXY3K/5hj/S7fvcvQzUzj4kq4CXgT+b0Tsn7pMwNeBfcCbwF9ExNGUA7W15Rh6R75Zucbfs/28LDLD/xxwHLhuxmX3AreN3j4MPD7612qWU+wd+TzkGH+HPw8bimwkaRfwCeCba2xyH/BErHoeuF7SzkRjtBku3b4rm9ife+81jn2mcvu/yeV3dqiKzvC/BnwR2LLG5TcCr018vDL63OnJjSQdAA4AbNqweaGB2rtyudHkFJL1vHGLar+OLb+O2q+jipxm/Z7tt2du8CXtB85GxBFJe9fabMbnrrgFRMQysAywdWl73reQDOUQ+hwj30TQq44hpwPC+P8wh/A7+s0qMsO/B/ikpH3AJuA6Sd+JiE9NbLMC3DTx8S7gVLphWtuxzyX0OcS9jFnjbvsgkEP4Hf1mzQ1+RDwKPAowmuF/YSr2AE8D/0XS91l9sPZcRJzGkmgz9m2GvqtxLyqXg0Db4fcST3NKn4cv6WGAiDgIHGL1lMxXWT0t86Ekoxu4IYa+75GfZ/Lnbzr+OYTf0a+XItq5W7l1aXvcve2BVq67C9qKfdOhH3rgF9H0AaCt8Dv663vuzGNHImJ3ma8tdFqmNWsIsX/jFjn2Cxrvs6b2W1v38tp+vKrP/NQKmWnjl72pG7YDn854X9Y9629rmcfr+vVw8DPSx9g78vVqas2/zfA7+uk4+JloOvZDCP3Fmy81fp3X/GZj49c51sSs/9x7r3H0O8zBb1nfZvVthL6NsK9lrbE0eSCoO/xtzPYd/TQc/Bb1aVbfVOhzivsiZo277oNAE+F39LvFwR+IumLfROi7Gvl5Jn+uOuNfZ/gd/W5x8FvS5Oy+jtjXHfq+Rn4t0z9vHQeAusLf9BKPo1+eg9+CpmLftVn90CK/nvG+6Fr4235CNluf//CqYV2OfR1/9HPx5kt/fLMr1bl/6jhwN/U3Hf7jrHIc/AZ1NfZ1ht6Kq2Of1fF/6+jny0s6Dely7FNqI/K37vpt7ddxcuWG2q9jrI7lnjduUdIlnqaWd7yevxgHvwGOfTOhbyLsi153nQeC1OF39PvPwe+JlLHvSujbDHxR02Os4wCQMvypH9B19PPi4Nesidl9jrGvK/RdiPx66jwApA5/16Jv8/lB2xo59mncuuu3f3zrmzp+rlQP7qa8p9fEA7l+EHc+z/A7LLfYpwx9H+O+nsmfN9Ws/+LNlyrP9rs20/fSzvo8w69J3bONvsa+rzP5RaTcBylm+12b6dva5gZf0iZJP5f0kqRXJH1lxjZ7JZ2TdGz09uV6htsNXbprWfXGnGr5wKG/UurwV5HyfP26o9+l21/TiizpXAQ+EhEXJC0BP5P0bEQ8P7XdTyNif/oh2rRUN5gUsa+qych/bMeJ5N/z8Jnbk3/PaeN9VHWpJ8UST1d4aWe2ucGP1Vc5vzD6cGn01s4rn3dAV5Zy2o593aGvI+5Fr6eug0CK8FeNfqo1fZ+5045CD9pKugo4Avx74L9HxAszNrtL0kvAKeALEfHKjO9zADgAsGnD5tKDzpVjX0xdsW8q8vNMjyP1AeDWXb919AvwLP9KhYIfEX8APijpeuBJSe+PiF9ObHIUuGW07LMPeAq4bcb3WQaWAbYubfe9hBa0Gfs6Qp9L5NczOcZU8a862696zn5Xom+XW+gsnYj4PfBj4ONTnz8fERdG7x8CliRtSzXILujC7L4vsf/YjhN/fOua1OOuul+r/J/m8LrF8/gB3MsVOUvnhtHMHknXAh8FTkxts0OSRu/vGX3f19MPd5hyOJUth9h3NfKzpPxZqp7N03b0c/j9HooiSzo7gW+P1vE3AD+IiGckPQwQEQeB+4FHJL0NvAU8OHqwdxC6MIuocsNsO/Z9ifws458txVJP1bX9vvJa/rvUVpe3Lm2Pu7c90Mp1p1Zn8Nteyikbe4e+nBThLxv9Kg/kpljPr3Mtv0/Bf+7MY0ciYneZr/Vf2vZcF2Pfp6WbRaX42cvu/7aXdurUhXvhTfBz6VSU++y+rDZj37T7rztaaLv/ef6Omkfyro/tOFFptl92eafNP87yGTv1c/B7LPdZ17S6Y1807It+fV0Hgq5FP/ULqFh6Dn4Fnt1frsrsvq7YV438oteROv5VH9Tt2gO5dc7y/eCtg99bZWf3Tce+jtA3Efki150y/lVn+4vyLL+f/KBthrp2XnIusb//uqOtxn5a6vGU3V9dexbSrv3+d4mD30NNz+7bllvop3U1+mV/H3J+7GjoZ+t4Saekof/ijLU5u8858tPGY23yTB+zaZ7hG1BuNufYLy7FuIcwy/eyTj0c/MxU/UXP+e40DDv2YymWoIb6h2lWjYNvpbT1QGDXY59Smeg3OcvP1ZCXYx38Evr2C9PUDbrqrLRvse/bzzMt93ubQ+Tg90ifb2B9jWMbSztdOU3T6/jpOfgZ6coveJlgVJnd9zX2Y134+fq2rDNUDv7A+YZsNhwOvlkGqszyfcaOFeXgW+28nJOfrqzjW1oOvi3EobBF5HoiQd/OtCuqyIuYb5L0c0kvSXpF0ldmbCNJ35D0qqRfSPLfjzcs1xtWFUOb3Q/t57XmFXkunYvARyLigqQl4GeSno2I5ye2uRe4bfT2YeDx0b9mZpaJuTP8WHVh9OHS6G36Ca/vA54Ybfs8cL2knWmHamZmVRRaw5d0laRjwFngcES8MLXJjcBrEx+vjD43/X0OSHpR0ouX3nmr7JjNzKyEQsGPiD9ExAeBXcAeSe+f2mTWAvIVL3sTEcsRsTsidm/ccO3iozUzs9IWOksnIn4P/Bj4+NRFK8BNEx/vAk5VGpmZmSVV5CydGyRdP3r/WuCjwPSJ1U8DnxmdrXMncC4iTicfra2pj68jOrQXCxnaz2vNK3KWzk7g25KuYvUA8YOIeEbSwwARcRA4BOwDXgXeBB6qabzWspMrN/hcfOu8jSdW2h5CK+YGPyJ+AXxoxucPTrwfwGfTDs3Muq6P9zy7zH9pa7U7fOb20l87lGWOKj9nmf17cuWG0tdn3eXgD9w1v9nY9hDMrCEOfka2/upi20MopMzs0LP8tXXh5/PEoB8c/B7p83ppF6JYRtWfq8/LOV2ZAHWJg1/CUB/hr6rKLB/6F/2+/TzT+jwB6SoH30rdXe/KLLHPqh5Ai+rbcs6QJ2wOfs/kPqtKMcvv+sw4xc9Qdj/6QD1sDn5m2lq3bHKWn2Jm2tXotznuMv9fZWf3VSceXr+vR5G/tLUZNp5YGeyr5uRiHM8uvHBIytA3tZRj/eMZfg+VnV11bZY/lvtsP4fYN7mUk/Oy4pDX78HBz1LX7s7mEv3cwp96TE3P7Nt6sLZrv/9d4iWdntry6yj1OrfX/GYjF2++tPDXlX1StcNnbudjO6affLW8ycC2sdRT10GnSuy7tHZv9XLwK6hzHX/rry5y7r3X1PK95ykb/bLGMUsZfmgu/nXfs2g69m2qc3Y/9OUccPB7rewsv6yqT52cerY/aa0oL3ogaHLZqOoSTtnYe3bfXw5+xro4y885+rPktu4/1rXYWzf4QduKcr+bWGXWVfbGX3UZ4fCZ2wd96mEXY59idu/lnPo5+JlLcSPoYvRheOFP8fO2sWafe+ztXQ5+An2fPbQZfeh/+FP9fFX2d5+Xcvp++1xEkRcxv0nSP0o6LukVSZ+bsc1eSeckHRu9fbme4Q5T27P8Kk6u3JA0/H2S8kDWVuw9u++WIg/avg38ZUQclbQFOCLpcET889R2P42I/emH2A11P9VCigdwq5y1U/VUzVQvfj4ZyCYf3E0l9UGr6sG07djXzbP7yxV5EfPTwOnR+29IOg7cCEwH3zqgavSB0uEfxylF+KE78a/rnkkfYu/ZfbMWOi1T0q3Ah4AXZlx8l6SXgFPAFyLilRlffwA4ALBpw+ZFx5q9Lszyofr5+bnM9idNR7XtA0Cdy08plsiGEHvP7q9UOPiSNgM/BD4fEeenLj4K3BIRFyTtA54Cbpv+HhGxDCwDbF3anv/9wRIc/WJSz/anNXkAaOqxhbZDb91XKPiSlliN/Xcj4kfTl08eACLikKTHJG2LiN+lG6ql1nb0of7wj3X9Ad9cYu/ZfbcVOUtHwLeA4xHx1TW22THaDkl7Rt/39ZQD7ZK6f9lS3liq3oBTzRhTns3TF+N94tgvxrFfW5EZ/j3Ap4GXJR0bfe5LwM0AEXEQuB94RNLbwFvAgxHRyyWbXKR82oUUM30o/2DupKZm/DlLeeDLKfTgB2nbpra6vHVpe9y97YFWrrspTbwiVsrn2knxRGt1PMvmUOKf+h7OEGM/hNn9c2ceOxIRu8t8rZ88rUZNvAxiTjN9SDvbH+vzrL+OZaxUy2yOff84+D2QW/Sh3vBDt+Nf52MVQ429FePg16ypFztPHX1Is8RT14upTEcz5wNAEw9Gpzzdsoux9+y+GAe/AV2MPuQ925+2VlSbPBC0cZZRrqEHxz5HDn5Dhh59uDxOTb2EYpEIFzko5HbKaOo/oHLsh8HBb1CXow9plnjGmpj1F5VbzNdS11/JOvbD4efDb1hTv6R13OjqeHbEa36z0X/uP0dd+2jLr6OzsbdyPMNvQVdn+lDPbB/aWe7JWZ0Hwbqe1rjJ2Ht2X46D33PjG2FXwg/DjX/d93T6EHpw7Ktw8FvS1Cx/rI7ZPtQbfrgygn06ADS5lOXYGzj4repL9KH+8I91+QDQxmMVdb4qlWPfPQ5+y9qIPqRf4hlrKvxja0W0zQNB2w9C1/3Sg208MOvYp+HgZ6Dp6EO9s31oPvzTikZ3kQND2yGfp4nXmHXsu83Bz8T4l7pPs31oP/zz5B7xIvoaenDsU3PwM9PH2T5cHqVc498lTUR+zLHvDwc/Q21FH+qd7Y85/uU0GXlw6PvIwc9UG9GHZsMPjv88TUce2v1rWce+Xg5+xtqKPjQffrgybkM8ALQR+LG2nxbBsa/f3OBLugl4AtgBvAMsR8TXp7YR8HVgH/Am8BcRcTT9cIenjQdzJ7UR/rEhHADaDPxY26EHx74pRWb4bwN/GRFHJW0Bjkg6HBH/PLHNvcBto7cPA4+P/rVE2pztQ7vhH5sVxy4dBHKI+ySHfnjmBj8iTgOnR++/Iek4cCMwGfz7gCdi9RXRn5d0vaSdo6+1RNqOPuQR/knzItr0ASG3qM+SQ+jBsW/DQmv4km4FPgS8MHXRjcBrEx+vjD53WfAlHQAOAGzasHmxkRrQ/hLP2GQ0con/LF0IcFMceiscfEmbgR8Cn4+I89MXz/iSK25pEbEMLANsXdruW2IFOcz2x3Kb9du7con8mGPfrkLBl7TEauy/GxE/mrHJCnDTxMe7gFPVh2fryWW2P9aVWX/f5RZ5cOhzMfcVr0Zn4HwLOB4RX11js6eBz2jVncA5r983Z+OJlexuUFt/dTHL8PRZrvs8t9/NISsyw78H+DTwsqRjo899CbgZICIOAodYPSXzVVZPy3wo/VBtnpyWecamA+SZfzo5xn2SQ5+fImfp/IzZa/ST2wTw2VSDsvJyW+aZ5mWfanKPPDj0OfNf2vZU7uEHz/6L6ELgJzn2eXPwe64L4R+bFbehHQS6Fvgxh74bHPyByHF9v4i+HgS6GvZpDn23OPgD0qXZ/nrWi2VOB4O+RH0Wh76bHPwBmryxdj3+0/oc2Rw49N3m4A9cX2b9Vh9Hvj8cfAMcfruSQ98/Dr5dps/LPVaMQ99fDr6tybP+4XDkh8HBt7k86+8nR354HHxbiOPfbY78sDn4Vprj3w2OvI05+JbEdFR8AGiXI2+zOPhWC8/+m+XAWxEOvtVuVox8EKjGgbcyHHxrhZeAinPcLRUH37LgewGrHHerk4Nv2ZoXv64eEBx1a8vc4Ev6W2A/cDYi3j/j8r3A/wL+dfSpH0XEf0s5SLNZFglnnQcHB9y6osgM/++AvwGeWGebn0bE/iQjMquBo2wGG+ZtEBE/Af6tgbGYmVmN5ga/oLskvSTpWUnvS/Q9zcwsoRQP2h4FbomIC5L2AU8Bt83aUNIB4ADApg2bE1y1mZkVVXmGHxHnI+LC6P1DwJKkbWtsuxwRuyNi98YN11a9ajMzW0Dl4EvaIUmj9/eMvufrVb+vmZmlVeS0zL8H9gLbJK0AfwUsAUTEQeB+4BFJbwNvAQ9GRNQ2YjMzK2Vu8CPiz+dc/jesnrZpZmYZS3WWjpmZZc7BNzMbCAffzGwgHHwzs4Fw8M3MBsLBNzMbCAffzGwgHHwzs4Fw8M3MBsLBNzMbCAffzGwgHHwzs4Fw8M3MBsLBNzMbCAffzGwgHHwzs4Fw8M3MBsLBNzMbCAffzGwg5gZf0t9KOivpl2tcLknfkPSqpF9IuiP9MM3MrKoiM/y/Az6+zuX3AreN3g4Aj1cflpmZpTY3+BHxE+Df1tnkPuCJWPU8cL2knakGaGZmaVyd4HvcCLw28fHK6HOnpzeUdIDVewEAF58789jMZaLMbAN+1/YgCvA40+rCOLswRvA4U/sPZb8wRfA143Mxa8OIWAaWASS9GBG7E1x/rTzOtDzOdLowRvA4U5P0YtmvTXGWzgpw08THu4BTCb6vmZkllCL4TwOfGZ2tcydwLiKuWM4xM7N2zV3SkfT3wF5gm6QV4K+AJYCIOAgcAvYBrwJvAg8VvO7lEuNtg8eZlseZThfGCB5naqXHqYiZy+1mZtYz/ktbM7OBcPDNzAai1uB35WkZCoxzr6Rzko6N3r7c9BhH47hJ0j9KOi7pFUmfm7FNq/u04Bhb35+SNkn6uaSXRuP8yoxtWv/9LDjO1vfnxFiukvRPkp6ZcVnr+3NiLOuNM4v9KemkpJdHY7jiVMxS+zMiansD/hS4A/jlGpfvA55l9Vz+O4EX6hxPhXHuBZ5pY2xT49gJ3DF6fwvwf4D/mNM+LTjG1vfnaP9sHr2/BLwA3JnTvlxgnK3vz4mx/Ffge7PGk8P+LDjOLPYncBLYts7lC+/PWmf40ZGnZSgwzixExOmIODp6/w3gOKt/1Typ1X1acIytG+2fC6MPl0Zv02cwtP77WXCcWZC0C/gE8M01Nml9f0KhcXbFwvuz7TX8tZ6WIUd3je5WPyvpfW0PRtKtwIdYnfFNymafrjNGyGB/ju7WHwPOAocjIst9WWCckMH+BL4GfBF4Z43Ls9ifzB8n5LE/A/gHSUe0+rQ00xben20Hv/DTMrTsKHBLRHwA+GvgqTYHI2kz8EPg8xFxfvriGV/S+D6dM8Ys9mdE/CEiPsjqX4fvkfT+qU2y2JcFxtn6/pS0HzgbEUfW22zG5xrdnwXH2fr+HLknIu5g9RmJPyvpT6cuX3h/th38TjwtQ0ScH9+tjohDwJKkbW2MRdISqyH9bkT8aMYmre/TeWPMaX+OxvB74Mdc+TTgre/LSWuNM5P9eQ/wSUknge8DH5H0naltctifc8eZyf4kIk6N/j0LPAnsmdpk4f3ZdvA78bQMknZI0uj9Pazut9dbGIeAbwHHI+Kra2zW6j4tMsYc9qekGyRdP3r/WuCjwImpzVr//Swyzhz2Z0Q8GhG7IuJW4EHgf0fEp6Y2a31/FhlnDvtT0nskbRm/D/wZMH0W4cL7M8WzZa5J9T0tQ9PjvB94RNLbwFvAgzF6mLxh9wCfBl4erekCfAm4eWKsbe/TImPMYX/uBL4t6SpWb9A/iIhnJD08Mc6292XRceawP2fKcH/OlOH+/BPgydFx52rgexHxXNX96adWMDMbiLaXdMzMrCEOvpnZQDj4ZmYD4eCbmQ2Eg29mNhAOvpnZQDj4ZmYD8f8BL5CqgL5MCLkAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] @@ -298,30 +314,152 @@ } ], "source": [ - "m = dist.Normal(torch.zeros(2), scale=torch.ones(2))\n", - "s = m.sample((1000,))\n", + "x1 = np.linspace(1, 5, num=50)\n", + "x2 = np.linspace(1, 5, num=50)\n", + "x1_s, x2_s = np.meshgrid(x1 ,x2)\n", "\n", - "plt.scatter(s[:, 0], s[:, 1], alpha=0.2)\n", - "plt.axis('scaled');" + "base = stats.multivariate_normal(mean=[3, 3], cov=0.5)\n", + "assert base.pdf([3, 3]) > base.pdf([1, 1])\n", + "\n", + "# field of all possible events X\n", + "x_field = np.concatenate([x1_s[..., None], x2_s[..., None]], axis=-1)\n", + "print(x_field.shape)\n", + "\n", + "plt.contourf(x1_s, x2_s, base.pdf(x_field))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Jacobian matrix.\n", + "\n", + "Assume the function $f: \\mathbb{R}^n \\mapsto \\mathbb{R}^m$. The jacobian is a 2D matrix, that stores on the partial derivatives of all the outputs and all the inputs. The width of the jacobian is equal to the functions input $n$. The height of the jacobian is equal to the functions output $m$.\n", + "\n", + "$$ \\mathbf{J} = \\begin{bmatrix}\n", + "\\frac{\\partial f_1}{\\partial x_1} & \\dots & \\frac{\\partial f_1}{\\partial x_n} \\\\[6pt]\n", + "\\vdots & \\ddots & \\vdots \\\\[6pt]\n", + "\\frac{\\partial f_m}{\\partial x_1} & \\dots & \\frac{\\partial f_m}{\\partial x_n} \\\\[6pt]\n", + "\\end{bmatrix}$$\n", + "\n", + "We will map 2D events to another domain. So we will have two inputs $n$ and two outputs $m$." ] }, { "cell_type": "code", - "execution_count": 23, + "execution_count": 275, + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/opt/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:7: UserWarning: To copy construct from a tensor, it is recommended to use sourceTensor.clone().detach() or sourceTensor.clone().detach().requires_grad_(True), rather than torch.tensor(sourceTensor).\n", + " import sys\n" + ] + } + ], + "source": [ + "def f(x1, x2):\n", + " return torch.exp(x1 / 3), x2**2\n", + "\n", + "def f_i(y1, y2):\n", + " return 3 * torch.log(y1), y2**0.5\n", + "\n", + "x_field = torch.tensor(x_field)\n", + "# Transform x events to y events\n", + "y_field = np.concatenate(f(x_field[..., 0, None], x_field[..., 1, None]), axis=-1)\n" + ] + }, + { + "cell_type": "code", + "execution_count": 276, "metadata": {}, "outputs": [], "source": [ - "def f(x):\n", - " torch.exp(x) + 2\n", - " \n" + "def create_det_jac(y_field):\n", + " # create for every y1, y2 combination the determinant of the jacobian f_i(y1, y2)\n", + " det_jac = np.zeros((y_field.shape[0], y_field.shape[1]))\n", + " for i in range(y_field.shape[0]):\n", + " for j in range(y_field.shape[1]):\n", + " y_field = torch.tensor(y_field)\n", + " y_field.requires_grad_(True);\n", + " fiy = torch.cat(f_i(y_field[..., 0, None], y_field[..., 1, None]), dim=-1)\n", + " fiy[i, j].sum().backward()\n", + " \n", + " # Ouputs of the partial derivatives are independent. I.e. f1 is dependent of y1 and not y2\n", + " # and vice versa f2 is dependent of y2 and not y1\n", + " # therefore the multiplication w/ 0 \n", + " row1 = y_field.grad[i, j].data.numpy() * np.array([1., 0.])\n", + " row2 = y_field.grad[i, j].data.numpy() * np.array([0., 1.])\n", + " \n", + " det = np.linalg.det(np.array([row1, row2]))\n", + " det_jac[i, j] = det\n", + " return det_jac\n", + "\n" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 277, "metadata": {}, - "outputs": [], - "source": [] + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/opt/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:5: UserWarning: To copy construct from a tensor, it is recommended to use sourceTensor.clone().detach() or sourceTensor.clone().detach().requires_grad_(True), rather than torch.tensor(sourceTensor).\n", + " \"\"\"\n" + ] + }, + { + "data": { + "text/plain": [ + "" + ] + }, + "execution_count": 277, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD8CAYAAABuHP8oAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAWAElEQVR4nO3df6ik1X3H8c/XuGuCGmNYfyyuuiCh2zQQYxexEYJtakjTUg0kJUJT2wY2LSkYEiipfzS2pZA/0rS0kBZTQyxNLKFqEoKRLNZWAm3Kam007NIkZbUbV1ebVlcqLuq3f9xn4uzs/Hh+nPOcH8/7BYv3zsy9c+7xzvueOfPMjLm7AADlOS31AAAA/RBwACgUAQeAQhFwACgUAQeAQhFwACjUxoCb2cVmdr+ZHTSz75rZTc3pt5jZD83s4ebfe+IPFwAwY5uOAzeznZJ2uvtDZna2pAclXS/pVyQ97+6fjj9MAMCi0zddwN2PSjrafHzczA5Kuij2wAAA621cgZ90YbPdkh6Q9BZJH5P065Kek3RA0sfd/X+WfM0+Sfsk6TU6/afPPP3coWMGgEl57qWnn3H38xZPbx1wMztL0j9J+mN3v8vMLpD0jCSX9Efa2mb5zXXf45xt5/vbd7y/8+ABYMruffKzD7r73sXTWx2FYmbbJN0p6YvufpckuftT7v6yu78i6XOSrgw5YADAem2OQjFJt0k66O6fmTt959zF3ivp0fDDAwCssvFBTElXS/qgpEfM7OHmtJsl3WBml2trC+WwpA9HGSEAYKk2R6F8S5ItOeue8MMBALTFMzEBoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKRcABoFAEHAAKtTHgZnaxmd1vZgfN7LtmdlNz+hvNbL+Zfa/577nxhwsAmGmzAn9J0sfd/SclXSXpI2b2ZkmfkHSfu79J0n3N5wCAkWwMuLsfdfeHmo+PSzoo6SJJ10m6vbnY7ZKujzVIAMCpOu2Bm9luSW+T9G1JF7j7UWkr8pLODz04AMBqrQNuZmdJulPSR939uQ5ft8/MDpjZgROvvNBnjACAJVoF3My2aSveX3T3u5qTnzKznc35OyUdW/a17n6ru+91973bT3tdiDEDANTuKBSTdJukg+7+mbmzvibpxubjGyV9NfzwAACrnN7iMldL+qCkR8zs4ea0myV9StKXzexDkh6X9P44QwQALLMx4O7+LUm24ux3hh0OAKAtnokJAIUi4ABQKAIOAIUi4ABQKAIOAIUi4ABQKAIOAIUi4ABQKAIOAIVq81R6FOLEnl2jXM/2Q0dGuR4A6xHwQo0V6y7XTdiBcRHwAqSMdReEHRgXAc9UKdFuY9nPQtSB4Qh4RmqK9ibzPysxB/oh4IlNKdqrLM4BQQfaIeCJEO7VWJ0D7RDwkRHubmbzRciBUxHwkRDuYViVA6ci4JER7vCIObCFgEdCuMdBzDFlvBZKBMQ7jRN7djH3mBQCHhgBSY+QYyrYQgmEYOSHI1hQOwIeAPHOG/vkqBUBH4Bwl4dVOWrCHnhPxLts7JOjBgS8B2749SDkKBkB74gbe50IOUpEwDvgBl4/Qo6S8CBmSyXfqJ+97IxeX3fOD14MPJJy8GAnSkDAW8g93n0DHeL71h75E3t2EXFki4BvkGO8YwW7j2VjqS3qrMaRKwK+Ri7xzinYbSyOt5agE3LkhoCvkDrepUV7ndqCzrYKckHAl0gZ75rCvcr8z1hqzFmNIwcEfEGKeE8h2quUHnNCjpQI+Jyx4z3lcC9TcszZVkEKG5/IY2afN7NjZvbo3Gm3mNkPzezh5t974g4zPuKdl2cvO6O4OeJJQBhbm2difkHSu5ec/qfufnnz756ww6pXiWFKaTZfJc0ZEcdYNm6huPsDZrY7/lDSGeMGV1KAcjWbwxK2V9gbxxiGvBbK75jZd5otlnNXXcjM9pnZATM7cOKVFwZcXRzEuzwlrchZjSOmvgH/S0mXSbpc0lFJf7Lqgu5+q7vvdfe92097Xc+riyP2jauk0JSolO0V9sYRS6+Au/tT7v6yu78i6XOSrgw7rPLlHpXalBJyIKRehxGa2U53P9p8+l5Jj667fI5i3phyCMnxSy3o9zv7MQ/6/WLJfZ+cww0R0saAm9kdkq6RtMPMjkj6pKRrzOxySS7psKQPRxxjcLHinSrcoWPd5TpyDfuzl52RdcQlHuDEcG2OQrlhycm3RRhL0caM9xjBbmtxLDkFndU4aje5Z2LGWH2PEe+cor1OjkHPOeREHENMKuClxbuUaK8z/zOkjnmuIWdLBX1NKuChxYp3DeFeJpfVea7746zG0dVk3tQ49Oo7RryPX2rVxnuZlD9vrocdcqghumAF3kPoG/6Uor1Mym2WHLdVWImjrUmswEOuakLGe2or7jZSzUluq3GevYk2JhHwUELHG6ulCHmO2ypEHOtUH/DcbgCsurtJFfKc5PY7jHywB95SiBt1qnC/eMmJXl93xuPbA4+kv9ncjbVHntuRKuyLY5mqAx5q5VJKvPuGus/3SxX3MUOe2wOcRByLqg54CEPjHTPcoYM95LrHDvrYISfiyFG1AQ+x+s4x3imjvU6qoI8V8pxW40QcM9U/iJlKyHi/eMmJH/8rxdhjHuvxhVwe4OSBTUiVBjz16jtUTEqL9ipjxXysI1aIOHJRZcCHSh3vWsK9zFghj42IIwfVBTzlL/TQcNQc7kWxf9YxVuNEHKlVF/Ch+t4oQ8R7imJvr4wR8RxCTsSniYAHMCQSU1p1bxJrLqayGifi01PVYYRDf4HHvhGOGe7du57udPnDR86LNJLNXrzkRJTDEI9falEPN8zheHEOMZyWqgI+xNhbJ7Hi3TXUXb/PWGGfzU/okBNx1KSagKe4+9gn3qHDHSrYQ64vZtRjhDz2k3+IOMbCHrj6rb5Tx3v3rqdHj/cqs7HEHE+M/fGY++LsiWMM1azAx5Qq3rkEe535McZYmYfeH4+5pZLDShx1q2IFPmSlMcZKaWi8c1ptdxFr3KzE22MVXrcqAj6mrjf2IbEpNdyLYvwcobdUiDhKNOmAd71hjR3v2sQKeSg1v1MSEa9T8QGv7RezllX3OqF/xpCr8VgRT70Kl+q7raCCgI9ljNV37eFeFCPkIdQccdRlskehxLwx5Rrvay881Opy+5/cE3kkJ5v97CGOWgl1lEqso1NSH5nC8eF1mWzAu+iyIssh3m1D3fXrY4d9966niTjQQdEB77unl9Nd2VDxHhrtPtcRI+ihVuNEfDVW4fUoOuBjiLn6HhrvMaLd9vpDxzzEajz3iKdExOvAg5iBjBnvay88lDzei2KMKcS9k5wf2MzpniDKNLmAd7nRxDoaYWi8cxY65EQ8Hg4rLF+xWyg5/fJ1CUTfIOUe7kWz8YbYWgmxLx7rNcZLx1ZK2Sa3Ak9pKvGeF3JFPnQ1HmIlXtsqHGXbGHAz+7yZHTOzR+dOe6OZ7Tez7zX/PTfuMMOIsX0S8111ctzr7ouIr8dWCvposwL/gqR3L5z2CUn3ufubJN3XfI41uoanlnDPC/UHaWrPWAVW2Rhwd39A0o8WTr5O0u3Nx7dLuj7wuNaKvWIIvfom3idLHXFW4adiFV6mvnvgF7j7UUlq/nv+qgua2T4zO2BmB0688kLPq5uO2uM9k/rnzDXiQBfRj0Jx91sl3SpJ52w7P9mzIUKvbmKtvkN73+sfGvT1f//cFYFGcqprLzw06CiVoU/2yfHIFJ6hiS76BvwpM9vp7kfNbKekYyEHlVLKVVWoVenQaK/6XjFinjriQ9X4LE2Uo+8Wytck3dh8fKOkr4YZTl26rL5DxPt9r38oaLzH+v5Df/bU++GhsReOttocRniHpH+W9BNmdsTMPiTpU5KuNbPvSbq2+XwyQt/oQ8V7LDFCnnpPfAj2wpHKxi0Ud79hxVnvDDyWqoy19z1muFddd6itlSHbKUO2UtgLR6mKeyZmn7t4be+SplhJDVl5poz3vJDjGDIfKR8wrmkVzjZKOYoLeGq57JnmEu+Z3MbTVS7/X4EuCHgEbVeCfVebucYy1LhKXYWHxmukYBMCXphc4z2T+/jWGboKr2kbBWUg4AUpJY4hxlnyUSk1YB+8DAQ8sNjbJ6VI+cempm0UYB0C3gEPdI2rxD9yobdR2AfHOgS8kfv+ZSnbJ/NKHDN/pFESAo4qsY2CKag+4DneBS1xa6CvoavwKc0V0FX1AQeAWhHwApS4lwwgPgIOAIUi4ABQKAIOAIUi4AWI+b6UAMpVfcBzfFH8Ie8BWZqhf3ymNFdAV9UHHNOU8o2OgbEQcAAoFAFvnP2Yb7xMyvdNLHEfvMTtk6H/j9v8HnWR4xYg8kHAA2t7171PnEqKeMqxsn2CqSDgCC5EvHnwMq3th46kHgJaIOCFyX0Vnvv41slt+wTYhIB31OZGHnMbRdqKZI6hDDWmIavvmrZP2P/GJpMIeNsbQooV1JBY5RLxkH9QUsWb1fer2D4pR3EBL+WXa6yVYOqIh7z+VPveKY8uWoXVN9o4PfUASnTG49uDvvXW/if3DHrjgllEx3rZ2Rh/NIbGO+XWSU2rb5SFgC84+zEP9v6Yh4+c1/qtvYZGXDo5rKFjHnOlnzLerL5PVso9XGyZTMDP+cGLQd9eLfQqXAoT8ZllwW0b9TG3ZUqPN6tvpDSZgKfSZRUuhY34otT75fNC7HfXGG9W3+iiuAcxpfi/aG1vmG0j0DU0tT+JJXW8c0W80dWkVuCht1G66LMSl+p6V/Ycwh1qz5vXPEEOilyBjyH0KlzqF5/9T+4pfkUe6meoNd6psfouFwEPIHbEpTK3VUL+8ak53mydoK9it1C2HzqiE3t2df66LtsoIQ8pnNd1O2VmPoY5bq3E+CMTYq+beC9HvMtXbMDH0jbiXQ8r7BvxmVxiHuueQe3hlog3hhsUcDM7LOm4pJclveTue0MMqq0xVuFd9Im4pEEhl06NaKygj7GNE+roEuK9GvGuR4gV+M+6+zMBvk+2umyl9HmCz9DV+KIS98tDHhZIvDEVk91C6boKHyPi0vDVeElCH8tNuDdj9V2XoUehuKRvmtmDZrYvxIC6GvIL2fUG1eWGfcbj23sF5fCR86p8ksq80D9j37leptZ4bz90hHhXaOgK/Gp3f8LMzpe038wOufsD8xdowr5Pkl572lkDry69rkem9H3NlPnA1bAqj/FHKeQLUcU8tjuHeKNO5h7mF9fMbpH0vLt/etVlztl2vr99x/uDXN+iPg9mzvR5QLPP4YUhXvyqlJjHvBcR+hUEa111S8S7Fvc++dkHlx0k0nsFbmZnSjrN3Y83H79L0h8OGGMyfY5K6XOMeIhXMFwMYw5BH2PLJ8bLvta86paI9xQM2UK5QNLdZjb7Pl9y93uDjKqHvocUzowZcSnMalxaHc/QYU+1L0+4uyPc09E74O7+n5LeGnAsg5UScSl8yBeV/EBorDdZiP0aJsQbY5vsYYSr9I241G9ffD5WsWKeu5jvijPGC08RbqRSXcCHrsKl/s/UHPraKVOJ+RhvYzaVcEvEe8qqC7iUPuJSv9X4vMXIlRz0sd53cqyXeSXcyEWVAZfSRlwKF/KZEoKe4g2CpxZtiXDjVdUGPJShL3wVOuQzm2IZI/C5vIP7mG+oQLiRs6oDHmIVLr16Iw4Rcil8zJfJJbahjP0uOIQbJag64FK4iEvhXoZ27JiXKMXbluUUbYlwY7PqAy6Fj7g0bDU+j5infY/J3KItEW60N4mAS2EjLoUPuXRqyGoMeg5vCJxjtCXCje4mE3Dp1RtI7iGfWRW73MOeQ6Tn5RpsiWhjmEkFfCb0alw6ORIxYj6vTSBDRz63KG+Sc7Qlwo0wJhlwKU7EZ2KuytsqLbhD5R5siWgjvMkGXIqzpTJvzFX5lJQQ6xmijZgmHfCZ2CGXTo0OQW+npFjPEG2MhYDPibmtsoign6zEUM8QbKRCwBeMsRpfZlXAagl7yYFeRLCRCwK+wvyNdOyYz2sTvpSRrynMqxBs5IqAt5BqVd7WFCI6JoKNUhDwDnJZlSMcYo2SEfCeiHl5iDVqQ8ADWAwDQU+LUGMqCHgEBD0+Ig0Q8FGsig1hX41AA5sR8ITWRarmuBNnIAwCnqkhkYsZf+IL5IOAV4jIAtNwWuoBAAD6IeAAUCgCDgCFIuAAUCgCDgCFMvfx3jvRzJ6W9NiGi+2Q9MwIwwmBscbBWONgrOGNNc5L3f28xRNHDXgbZnbA3femHkcbjDUOxhoHYw0v9TjZQgGAQhFwAChUjgG/NfUAOmCscTDWOBhreEnHmd0eOACgnRxX4ACAFgg4ABQqScDN7PNmdszMHl1x/jVm9qyZPdz8+/2xxzg3lovN7H4zO2hm3zWzm5Zcxszsz83s+2b2HTO7IuOxJp9bM3utmf2rmf17M84/WHKZXOa0zViTz+nCeF5jZv9mZl9fcl4W8zo3nnVjzWZezeywmT3SjOPAkvPTzKu7j/5P0jskXSHp0RXnXyPp6ynGtmQsOyVd0Xx8tqT/kPTmhcu8R9I3JJmkqyR9O+OxJp/bZp7Oaj7eJunbkq7KdE7bjDX5nC6M52OSvrRsTLnMa8uxZjOvkg5L2rHm/CTzmmQF7u4PSPpRiuvuyt2PuvtDzcfHJR2UdNHCxa6T9De+5V8kvcHMdo481LZjTa6Zp+ebT7c1/xYfTc9lTtuMNRtmtkvSL0r66xUXyWJepVZjLUmSec15D/xnmrut3zCzn0o9GEkys92S3qatVdi8iyT919znR5Q4nGvGKmUwt81d54clHZO0392zndMWY5UymNPGn0n6XUmvrDg/m3nV5rFK+cyrS/qmmT1oZvuWnJ9kXnMN+EPaeu7/WyX9haSvJB6PzOwsSXdK+qi7P7d49pIvSbZK2zDWLObW3V9298sl7ZJ0pZm9ZeEi2cxpi7FmMadm9kuSjrn7g+sutuS00ee15VizmNfG1e5+haRfkPQRM3vHwvlJ5jXLgLv7c7O7re5+j6RtZrYj1XjMbJu2gvhFd79ryUWOSLp47vNdkp4YY2yLNo01t7l19/+V9I+S3r1wVjZzOrNqrBnN6dWSftnMDkv6O0k/Z2Z/u3CZXOZ141gzmle5+xPNf49JulvSlQsXSTKvWQbczC40M2s+vlJb4/zvRGMxSbdJOujun1lxsa9J+rXmkeirJD3r7kdHG2SjzVhzmFszO8/M3tB8/DpJPy/p0MLFcpnTjWPNYU4lyd1/z913uftuSR+Q9A/u/qsLF8tiXtuMNZd5NbMzzezs2ceS3iVp8Qi6JPOa5E2NzewObT3CvMPMjkj6pLYeHJK7/5Wk90n6bTN7SdILkj7gzUO9CVwt6YOSHmn2QSXpZkmXSD8e7z3aehT6+5L+T9JvJBin1G6sOcztTkm3m9lrtHWj/LK7f93MfmtunLnMaZux5jCnK2U6r0tlOq8XSLq7+VtyuqQvufu9OcwrT6UHgEJluYUCANiMgANAoQg4ABSKgANAoQg4ABSKgANAoQg4ABTq/wE/9OE0GYs16wAAAABJRU5ErkJggg==\n", + "text/plain": [ + "
" + ] + }, + "metadata": { + "needs_background": "light" + }, + "output_type": "display_data" + } + ], + "source": [ + "px = base.pdf(x_field)\n", + "transformed = px * np.abs(create_det_jac(y_field))\n", + "plt.contourf(y_field[..., 0], y_field[..., 1], transformed)" + ] + }, + { + "cell_type": "code", + "execution_count": 278, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "0.9907110850291531" + ] + }, + "execution_count": 278, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "np.trapz(np.trapz(transformed, y_field[:, 0, 1], axis=0), y_field[0, :, 0])" + ] } ], "metadata": {