ishay2b/finding-points-for-given-distances.ipynb

## finding-points-for-given-distances.ipynb
{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "Finding points for given distances.ipynb",
      "provenance": [],
      "collapsed_sections": [],
      "authorship_tag": "ABX9TyNyWO9Nsapp3Vz+azPgTTm8",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/gist/ishay2b/7aba5d8c7832fe63eedcffab33dd2714/finding-points-for-given-distances.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "MR5udthOizcS",
        "outputId": "958fbd75-8168-4ad7-e706-a7b921843b63"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Success reaching total distance 0.002 after 303 iterations\n",
            "Possible Found positions:\n",
            "[[3.9977682  0.52691495]\n",
            " [5.2908134  3.2339554 ]\n",
            " [0.38747454 2.249567  ]]\n",
            "Accumlated error 0.0013997554779052734\n",
            "Actual distnaces      : [3.0000055 4.000219  5.0011754]\n",
            "Vs requested distnaces: [3. 4. 5.]\n"
          ]
        }
      ],
      "source": [
        "import torch\n",
        "from numpy import array as _A, dtype\n",
        "from torch.autograd import Variable\n",
        "import torch.optim as optim\n",
        "import numpy as np\n",
        "import matplotlib.pyplot as plt\n",
        "\n",
        "N = 3\n",
        "device = \"cpu\"\n",
        "\n",
        "\n",
        "def pairs(n):\n",
        "    ''' Create all valid non repeating pairs '''\n",
        "    return [[a, b] for a in range(n) for b in range(a, n) if a != b]\n",
        "\n",
        "\n",
        "all_pairs = _A(pairs(N))\n",
        "src_index = all_pairs[:, 0]  # Index of the points p1\n",
        "trg_index = all_pairs[:, 1]  # Index of the points p2\n",
        "\n",
        "#The input distances vector, assume ordered same as the pairs function above.\n",
        "all_distances = torch.tensor([3.0, 4.0, 5.0], requires_grad=False, device=device)\n",
        "\n",
        "#Random starting points scaled to mean distance as nicer starting points.\n",
        "start_points = np.random.random((3, 2)) * all_distances.mean().detach().cpu().numpy()\n",
        "\n",
        "# The coordintes we are trying to find.\n",
        "positions = torch.tensor(start_points, requires_grad=True,\n",
        "                         device=device, dtype=torch.float)  \n",
        "\n",
        "optimizer = optim.Adam([positions], lr=0.01)\n",
        "\n",
        "loss_graph = []  # For visual\n",
        "threshold = 0.002  # Min accumulated error for early stopping.\n",
        "\n",
        "for i in range(1000):\n",
        "    optimizer.zero_grad()\n",
        "    p1 = positions[src_index]  # Vector of start points.\n",
        "    p2 = positions[trg_index]  # Vector of end points.\n",
        "    euc = ((p1 - p2)**2).sum(axis=1)**0.5  # Vector of euclidian distances.\n",
        "    loss = (euc - all_distances).abs().sum()  # Total loss error to minimize.\n",
        "    loss_graph.append(float(loss.detach().cpu().numpy()))\n",
        "    if loss_graph[-1] < threshold:\n",
        "      print(f\"Success reaching total distance {threshold} after {i+1} iterations\")\n",
        "      break\n",
        "    loss.backward()  # Now Optimize!\n",
        "    optimizer.step()\n",
        "\n",
        "\n",
        "print(\"Possible Found positions:\")\n",
        "print(positions.detach().cpu().numpy())\n",
        "print(f\"Accumlated error {loss_graph[-1]}\")\n",
        "print(\"Actual distnaces      :\", euc.detach().cpu().numpy())\n",
        "print(\"Vs requested distnaces:\", all_distances.detach().cpu().numpy())\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "plt.plot(loss_graph, label=\"loss_graph\")\n",
        "plt.legend();"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 265
        },
        "id": "1pm-HyUci7I9",
        "outputId": "b9031b33-6cac-4730-ebca-22e4a2a78ed1"
      },
      "execution_count": 3,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAD4CAYAAADFAawfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3dd3hUZcL+8e8zk0kmCUkIEHogiSBVpASkgyh2xcIuNlSau7bFuuVd9/fuvuv77q51sRdAQFGx4LriuisivUmIQdrSCYSaUAIJJKQ8vz8ysIgBJpDJmZncn+viYjIt95ODtyfPnHMeY61FRESCl8vpACIicmYqahGRIKeiFhEJcipqEZEgp6IWEQlyEYF40wYNGtiUlJRAvLWISFhavnx5nrU2qbLHAlLUKSkpZGRkBOKtRUTCkjEm+3SPaepDRCTIqahFRIKcilpEJMj5NUdtjKkLjAc6AhYYaa1dHMhgIuKMkpIScnJyKCoqcjpKWPJ6vTRv3hyPx+P3a/z9MHEc8E9r7VBjTCQQcy4BRST45eTkEBcXR0pKCsYYp+OEFWst+/btIycnh9TUVL9fd9apD2NMAtAfmOD7RsestQfPOamIBLWioiLq16+vkg4AYwz169ev8m8r/sxRpwK5wNvGmO+MMeONMbGVBLjXGJNhjMnIzc2tUggRCS4q6cA5l5+tP0UdAXQFXrPWdgEKgV+f+iRr7ZvW2nRrbXpSUqXHbJ9RUUkZb87bxJLN+6r8WhGRcOZPUecAOdbapb6vP6aiuKs3iDFMWLCFl7/ZWN1vLSIS0s5a1Nba3cB2Y0wb312XAWuqO0hkhIu7e6ewYGMea3Yequ63F5EQUqdOHacjnLOUlBTy8vKq9T39PY76IWCqMeZ7oDPwf9WawueOHi2JiXQzfsHmQLy9iIjfSktLnY5wgl+H51lrs4D0AGchIcbDT9OTeXdJNr+8si2NE7yB/pYicgZ/+Hx1tf+G275pPP99fQe/nmut5Ze//CVffvklxhiefPJJhg0bxq5duxg2bBiHDh2itLSU1157jd69ezNq1CgyMjIwxjBy5EgeeeSRSt932bJljBo1CpfLxeDBg/nyyy9ZtWoVkyZNYvr06RQUFFBWVsYXX3zBkCFDOHDgACUlJTz11FMMGTKErVu3ctVVV9GtWzcyMzPp0KEDU6ZMISam4sjll156ic8//5ySkhI++ugj2rZte14/s6A7M3FU31TKrWXSoq1ORxERh02fPp2srCxWrFjB119/zRNPPMGuXbt47733uPLKK0881rlzZ7KystixYwerVq1i5cqVjBgx4rTvO2LECN544w2ysrJwu90/eCwzM5OPP/6YuXPn4vV6+fTTT8nMzGT27Nk89thjHF9ndt26ddx///2sXbuW+Ph4Xn311RPv0aBBAzIzM7nvvvt49tlnz/vnEJCr552P5HoxXNWxMe8tzeahQa2IjQq6iCK1hr97voGyYMECbrvtNtxuN40aNWLAgAEsW7aM7t27M3LkSEpKSrjxxhvp3LkzaWlpbN68mYceeohrr72WK664otL3PHjwIIcPH6ZXr14A3H777cyYMePE44MHD6ZevXpAxR79f/3XfzFv3jxcLhc7duxgz549ACQnJ9OnTx8A7rzzTl588UUef/xxAG6++WYAunXrxvTp08/75xB0e9QAo/ulcaiolA8ztjsdRUSCUP/+/Zk3bx7NmjXjnnvuYcqUKSQmJrJixQoGDhzI66+/zujRo8/pvWNj/3OayNSpU8nNzWX58uVkZWXRqFGjEyernHo89MlfR0VFAeB2u6tlrjsoi7pri0TSWyYyYcEWSsvKnY4jIg7p168f06ZNo6ysjNzcXObNm0ePHj3Izs6mUaNGjBkzhtGjR5OZmUleXh7l5eXccsstPPXUU2RmZlb6nnXr1iUuLo6lSyuOOP7ggw9O+/3z8/Np2LAhHo+H2bNnk539n0tGb9u2jcWLKy559N5779G3b99qHPkPBe28wuh+afz83eX8a/Ueru3UxOk4IuKAm266icWLF3PxxRdjjOHpp5+mcePGTJ48mWeeeQaPx0OdOnWYMmUKO3bsYMSIEZSXV+zc/elPfzrt+06YMIExY8bgcrkYMGAACQkJlT7vjjvu4Prrr+eiiy4iPT39Bx8KtmnThldeeYWRI0fSvn177rvvvuod/EnM8Ynx6pSenm7Pd4WXsnLLoOfmkBgTyaf399YprSI1ZO3atbRr187pGAFVUFBw4ljtP//5z+zatYtx48b5/fqtW7dy3XXXsWrVqnP6/pX9jI0xy621lR5dF5RTHwBul2FU31Syth9kefYBp+OISBj54osv6Ny5Mx07dmT+/Pk8+eSTTkc6o6Cd+gAY2q05z89cz1vzN5OeUs/pOCISYh544AEWLlz4g/vGjh3LiBEjGDZs2Dm/b0pKyjnvTZ+LoC7qmMgI7rykJa/M2ciWvEJSG/zoon0iEgDW2rCYbnzllVecjvAj5zLdHLRTH8fd1bslHpeLiQu2OB1FpFbwer3s27fvnApFzuz4wgFeb9XOug7qPWqAhnFehnRuykfLt/Po4AtJjI10OpJIWGvevDk5OTnouvKBcXwprqoI+qIGGNM/jY+W5zB1aTYPDmrtdByRsObxeKq0TJQEXtBPfQBc2CiOARcmMWlRNkUlZU7HERGpUSFR1ABj+qWRV1DM37N2Oh1FRKRGhUxR92lVn7aN4xi/YLM+5BCRWiVkitoYw73901i/p4C56/Uhh4jUHiFT1ADXdWpKo/go3pqvFWBEpPYIqaKOjHBxT+9UFm7cx+qd+U7HERGpESFV1AC392hBTKSbCfN1AoyI1A4hV9QJMR6GdU/m7yt2sju/yOk4IiIBF3JFDTCyj9ZVFJHaIySLOrleDFd3bMLUpdkUFAfPku4iIoEQkkUNMLpfKoeLSvlwmdZVFJHwFrJF3cW3ruLEhVpXUUTCm19FbYzZaoxZaYzJMsac3xpb1WhM/zRyDhzlX6v3OB1FRCRgqrJHfam1tvPp1vRywuXtGpFSP4Y35+u0chEJXyE79QH/WVdxxfaDZGhdRREJU/4WtQW+MsYsN8bcW9kTjDH3GmMyjDEZNXnB8aHdkqkb4+GteTqtXETCk79F3dda2xW4GnjAGNP/1CdYa9+01qZba9OTkpKqNeSZREe6Gd6zJTPX7mFLXmGNfV8RkZriV1Fba3f4/t4LfAr0CGSoqhreq2JdxQkLtFctIuHnrEVtjIk1xsQdvw1cAdTcOul+aBjn5cYuTfl4eQ4HCo85HUdEpFr5s0fdCFhgjFkBfAt8Ya39Z2BjVd3ofmkUlZTz7pJsp6OIiFSrsy5ua63dDFxcA1nOy4WN4hjYJonJi7MZ0z8Nr8ftdCQRkWoR0ofnnUrrKopIOAqrou59QX3aNYnnLZ0AIyJhJKyK2hjDmH6pbNhbwBytqygiYSKsihoq1lVsHO9lvNZVFJEwEXZFHRnh4p4+KVpXUUTCRtgVNcBtPVoQG+lmvNZVFJEwEJZFnRDt4afdk/l8xU525R91Oo6IyHkJy6IGrasoIuEjbIs6uV4MV1/UhPeWbtO6iiIS0sK2qKHiBJjDRaVM07qKIhLCwrqoOyfXpXtKIhMXaF1FEQldYV3UULFXvePgUf65erfTUUREzknYF/Xl7RqR2iCWt+Zv0WnlIhKSwr6oXS7DSK2rKCIhLOyLGmBo1+Ykxnh4U+sqikgIqhVFHR3p5s6eLfla6yqKSAiqFUUNcFevFK2rKCIhqdYUdVJcFDd1acbHy3PYr3UVRSSE1JqiBhjdL1XrKopIyKlVRd3at67ilMVbKSopczqOiIhfalVRw/F1FY9pXUURCRm1rqh7X1Cfto3jGL9A6yqKSGiodUVdsa5iGuv3FDB/Q57TcUREzqrWFTXA9Rc3pWFcFG9pXUURCQG1sqgjI1zc3TuF+RvyWLf7sNNxRETOyO+iNsa4jTHfGWNmBDJQTbnjkhZEe9xarVxEgl5V9qjHAmsDFaSm1Y2JZGi35nyWtZO9h4ucjiMiclp+FbUxpjlwLTA+sHFq1si+qZSUl/PuYp0AIyLBy9896r8CvwROu0yKMeZeY0yGMSYjNze3WsIFWmqDWC5v14h3lmRz9JhOgBGR4HTWojbGXAfstdYuP9PzrLVvWmvTrbXpSUlJ1RYw0Mb0S+PAkRKmf5fjdBQRkUr5s0fdB7jBGLMV+AAYZIx5N6CpalD3lEQ6NU9gwvwtlJfrBBgRCT5nLWpr7W+stc2ttSnArcA31to7A56shhhjGNU3lc15hcxet9fpOCIiP1Irj6M+1TUXNaFpgpfx87c4HUVE5EeqVNTW2jnW2usCFcYpHreLe/qksHjzPlbtyHc6jojID2iP2mdY9xbERrqZsEB71SISXFTUPgnRHoZ1b8HnK3ayK/+o03FERE5QUZ9kRJ8Uyq1l8iKdACMiwUNFfZLkejFc3bEJ7y3NprC41Ok4IiKAivpHRvVL5VBRKR9lbHc6iogIoKL+ka4tEunaoi4TF26lTCfAiEgQUFFXYky/NLbtP8LMNbudjiIioqKuzBUdGpNcL1onwIhIUFBRV8LtMozsk0pG9gG+23bA6TgiUsupqE/jJ+nJxHkjGK8TYETEYSrq06gTFcHtPVrw5cpdbN9/xOk4IlKLqajP4J4+KbiMYdKirU5HEZFaTEV9Bk0Sorm2UxOmLdvOoaISp+OISC2loj6L0X3TKCgu5cNlOgFGRJyhoj6Li5oncElqPd5euJXSstMuGSkiEjAqaj+M7pfGjoNH+XKVToARkZqnovbDZW0bktoglvHzN2OtTisXkZqlovaDy2UY2TeVFTn5ZGTrBBgRqVkqaj8N7dqcujEexs/f7HQUEallVNR+io50c+clLflqzR625hU6HUdEahEVdRXc1aslHpeLtxfqtHIRqTkq6ipoGO/lhs5N+TAjh4NHjjkdR0RqCRV1FY3qm8rRkjLe+3ab01FEpJZQUVdRuybx9G3VgMmLtnKsVCfAiEjgnbWojTFeY8y3xpgVxpjVxpg/1ESwYDa6Xyp7DhUz4/udTkcRkVrAnz3qYmCQtfZioDNwlTGmZ2BjBbcBFybRumEdxs/fohNgRCTgzlrUtkKB70uP70+tbidjDKP7pbJm1yEWb9rndBwRCXN+zVEbY9zGmCxgLzDTWru0kufca4zJMMZk5ObmVnfOoDOkczPqx0ZqBRgRCTi/itpaW2at7Qw0B3oYYzpW8pw3rbXp1tr0pKSk6s4ZdLweN8N7teSbf+9l496Cs79AROQcVemoD2vtQWA2cFVg4oSW4T1bEhnh4q15Oq1cRALHn6M+kowxdX23o4HBwL8DHSwU1K8TxW3dk5mWsZ1HpmXp1HIRCYgIP57TBJhsjHFTUewfWmtnBDZW6Pjdde2Jj/bwxrzNfJa1g+s6NeWJK9uQXC/G6WgiEiZMIA4vS09PtxkZGdX+vsFs7+EiJi7YyqRFWyi38PP+aTwwqBVREW6no4lICDDGLLfWplf2mM5MrCYN47z8+uq2fPPYQK7u2JgXv9nI9S8tYGVOvtPRRCTEqairWdO60Yy7tQtvj+hO/tESbnx1Ic/PXK/1FkXknKmoA+TSNg356uEB3HBxU16ctYG7Jn5LXkGx07FEJASpqAMoIcbDC8M68/QtncjIPsD1Ly0gc5uW8hKRqlFR14Cfdk9m+n29iXAbhr2xmPHzN7Pz4FGnY4lIiNBRHzUo/0gJD76fyfwNeQD0vzCJX17Zho7NEhxOJiJO01EfQSIhxsOUkT2Y8VBfHh18Iat25DPklYW8NW+zrsInIqeloq5hxhg6NkvgF5e1ZvbjAxncrhH/+4+13D81k6KSMqfjiUgQUlE7KCHaw2t3duW317Tjn6t3M2ryMo4cK3U6logEGRW1w4wxjOmfxvM/vZjFm/Zx98RvOVxU4nQsEQkiKuogcVOX5rx0W1e+23aQO8YvZfv+I05HEpEgoaIOItd2asLrd3Zjc24hV4+bz8fLc5yOJCJBQEUdZC5v34gvx/ajfdN4Hv9oBb+ZvlKrnYvUcirqIJRcL4b3x/TkvoEX8P6327j9rSXsPVzkdCwRcYiKOki5XYZfXdWWl2/vwuqdh7jhpYWs233Y6Vgi4gAVdZC7rlNTPrmvN+XW8tM3FrM8W9cKEaltVNQhoH3TeD65rzeJMR7uHL+UuevDf5V3EfkPFXWISK4Xw0c/701qg1hGT17GV6t3Ox1JRGqIijqEJMVF8cHPetKhaQL3T81UWYvUEirqEBPv9TBlVA86NKso60kLt1Bergs6iYQzFXUIivd6eGdUDwZcmMTvP1/DtS8t4Lmv1rFhj44KEQlHKuoQFe/18NZd6Tz3k4uJdBtenbOJwS/M4xfvf8cOLUogEla0cECY2FdQzKRFW3lz3mYA/mdIB4Z1b+FwKhHxlxYOqAXq14nisSva8M3jA0lPSeRXn6zkiY9WcPSYrnEtEurOWtTGmGRjzGxjzBpjzGpjzNiaCCbnplndaKaMvIRfDGrFx5k53PTqQrL3FTodS0TOgz971KXAY9ba9kBP4AFjTPvAxpLz4XYZHr2iDW/f051d+UXc8PJCFvjWaRSR0HPWorbW7rLWZvpuHwbWAs0CHUzO38A2Dfn7g31oHO9l+MSljP3gOzK27tf6jCIhpkofJhpjUoB5QEdr7aFTHrsXuBegRYsW3bKzs6svpZyXwuJSXpy1gXeXZFN4rIxuLRP5/fUduKi5Vj8XCRZn+jDR76I2xtQB5gL/a62dfqbn6qiP4FRYXMr073bw0qwN7Cs8xph+aTx8eWu8HrfT0URqvfM+6sMY4wE+AaaeraQleMVGRTC8Z0tmPjqAoV2b8/rcTVw9bj5LN+9zOpqInIE/R30YYAKw1lr7fOAjSaAlRHv4y9BOTB19CaXl5Qx7cwlP/m0lhcVaAV0kGPmzR90HGA4MMsZk+f5cE+BcUgP6tGrAvx7uz+i+qUxduo3rXlrAqh35TscSkVPozEQBYOnmfTw8LYu8gmKeuLINI/ukEuHW+VAiNUVnJspZXZJWny/H9uOyto34v3/8m75/mc24rzeQf7TE6WgitZ6KWk6oGxPJa3d2Zfxd6VzYOI4Xvl7PoGfnMG3ZNl1KVcRBKmr5AWMMl7dvxJSRPZjxUF9SG8Tyq09WcuOrC1m9U/PXIk5QUctpdWyWwEc/78Vfh3VmV34RQ15eyLivN1BSVu50NJFaRUUtZ2SM4cYuzfjq4f5c26kJL3y9njvGL9U1r0VqkIpa/JIYG8m4W7vwwrCLWZmTz2XPzeGFmespKtFlVEUCTUUtVXJTl+Z8/dgALm/XiHGzNjD4hbnMXrfX6VgiYU1FLVXWrG40L9/elffH9CTS7WLE28u4f+pyducXOR1NJCypqOWc9bqgPl+O7c8TV7Zh1tq9XPbcHCYu2KIPG0WqmYpazktkhIsHLm3FzEcGkJ5Sj/+ZsYZef5rFhAVbKFVhi1QLFbVUixb1Y5g0ojuTR/agXZN4/jhjDcMnfKvpEJFqoKKWamOMYcCFSUwZ2YNnhnZi+bYD9H9mNn/4fDV7D6uwRc6VilqqnTGGn6QnM+vRAdzUuRlTFmfT/+mKa4ccK9V0iEhVqaglYJLrxfCXoZ34xnc43wtfr+f6lxbw3bYDTkcTCSkqagm4lvVjefn2rky4O538oyXc/Noi/jhjDUeOaaECEX+oqKXGXNauETMf7c8dl7RgwoItDH5+Hv9ctVuroouchYpaalSc18NTN17Exz/vRZw3gp+/u5wRk5aRva/Q6WgiQUtFLY5IT6nHjIf68uS17Vi2ZT+DX5jH81+t07qNIpVQUYtjItwuRvdLY9ZjA7mqQ2Ne/GYjlz47h799t0PTISInUVGL4xoneHnxti58cl9vmiR4eXhaFndN/JZNuQVORxMJCipqCRrdWiYy/f4+/HFIB7K2HeSKF+bxm+nfs+eQTpaR2k1FLUHF7TIM75XC7CcGMrxnSz5ensPlz83lnSXZWrdRai0VtQSlBnWi+P0NHZj5yAA6JSfwu7+t4pbXF2ndRqmVVNQS1FIaxPLuqEt4/qcXk73vCNe9tIBHp2WxU0uBSS1y1qI2xkw0xuw1xqyqiUAipzLGcHPX5sx+bCD39k9jxspdDH5+LuPnb6a4VEuBSfjzZ496EnBVgHOInFVCjIffXN2OWY9WXPv6qS/Wctlzc/l8xU4dzidh7axFba2dB+yvgSwifkmuF8PkkT14Z1QP4rweHnr/O37y+mKWZ+tiTxKeNEctIatf6yRmPNSXv9xyEVv3FXLLa4u4Y/wSMnV1Pgkz1VbUxph7jTEZxpiM3Nzc6npbkTNyuwzDurdg7hOX8ttr2rFudwG3vLaI3366kvyjJU7HE6kWxp+5PWNMCjDDWtvRnzdNT0+3GRkZ55dM5BwUFpfy/Mz1vL1wC/XrRDH2stYM6dyUOK/H6WgiZ2SMWW6tTa/sMU19SFiJjYrgd9e157MH+tIkwcuTf1vFpc/O4bMsXT9EQpc/h+e9DywG2hhjcowxowIfS+T8XNQ8gc8e6MMn9/Wmad1oxn6QxfAJ37J+z2Gno4lUmV9TH1WlqQ8JJmXllncWb+X5mespPFbGHZe04BeXtaZBnSino4mcoKkPqdXcLsM9fVKZ88Sl3NYjmalLtzHwmTm8/M0Gjh7TCTMS/FTUUmvUi43kqRsv4l8P96fXBfV59qv1mr+WkKCillqnVcM6vHVXOh/+rBcN46MY+0EWoydnsCtf1w+R4KSillqrR2o9Pr2/D7+7rj0LN+Ux4Jk5/PdnqzhQeMzpaCI/oKKWWs3tMozqm8rMRwZwS9dmvLt0GwOfncOkhVsoKSt3Op4IoKM+RH5g3e7D/HHGGhZszCO1QSwPXNqKIZ2b4nFrn0YCS0d9iPipTeM43hnVg/F3pRPtcfP4Ryu48q/zWLJ5n9PRpBZTUYucwhjD5e0b8cUv+vLm8G6UlJVz65tLeHRaFlvyCp2OJ7VQhNMBRIKVMYYrOjSmX+skxs3awNsLt/D3FTu5q1cKYy9rTUKMrh8iNUN71CJnER3p5tdXt2XBrwbxk/Rk3l60hYHPzuadJdmU6gNHqQEqahE/JcVF8aebL+KLh/rRpnEcv/vbKq4eN5/Z6/Y6HU3CnIpapIraN43n/TE9ecM3fz3i7WXcPfFbtmr+WgJERS1yDowxXNmhMV89MoAnr23H8uwDXP78XH776Uqd4SjVTsdRi1SDvYeKePGbDUxbth2D4Y6eLfjFoNYkxkY6HU1CxJmOo1ZRi1Sj7fuP8PI3G/lo+XbivB5+cVlrhvdsSWSEfnmVM9MJLyI1JLleDH8Z2okvx/anU/ME/jhjDVf+dR4zvt9Jebmu0CfnRkUtEgBtGscxZWQP3r6nO5FuFw++9x03vLKAuetzdUlVqTIVtUiAGGO4tG1D/jG2Hy8Mu5iDR0q4e+K33PbWEjK3HXA6noQQFbVIgLldhpu6NOebxwbyhxs6sHFvATe/uogxUzK0hqP4RR8mitSwwuJS3l64hTfmbqbgWCk3d2nOw5e3JrlejNPRxEE66kMkCB0oPMZrczcxadFWrLXccUlLHhzUSovu1lIqapEgtiv/KC/O2sCHGTlERbgY3S+NMf1SifPqok+1iYpaJARsyi3g+a/W88XKXSTGeLh/YCvu7NmS6Ei309GkBqioRULIypx8nv7Xv5m/IY8GdSIZ0y+NO3u2JDZKVyUOZypqkRC0bOt+Xpy1gfkb8kiM8TC6Xxp39WqpKZEwdd5FbYy5ChgHuIHx1to/n+n5KmqR6pO57QAvzdrA7HW51I3x8MDAVgzv1RKv54dTItZacguKOVJcRlFpGSWlluhIN80To3/0XAk+51XUxhg3sB4YDOQAy4DbrLVrTvcaFbVI9fs+5yDPfrWeeetzaZrg5e7eKfS/MIl563PJyD5AZvYB9hUe+9HrjIEm8V6aJUbTtUUiPdPq07tVfQ4UlrA8+wA5B46QV1BMVISb3YeKMEBaUh3aNYnjwkZxGAPFJeVERrg4XFRKfHQEDepEEeEyWAvl1lJYXMaRklIaxnnZX3iMtbsOERvlplndGOrFRrJ21yE27i3A7TJckFSHCxvXYU9+MbkFxXg9LkrLLAnRHprU9XLoaCml5eUcPFLCzoNH8XrcNK0bTZMELy5jOFxUwp5DxRhTMb6C4lIaxXmJ8rhYtSOfPYeKaZLgpXtqPQqLS9mUW4C1kJwYQ3K9aA4eKaGguBSvx43X48LrceNxu7DWUlJm2ZRbwO5DRTSrG03DuCg27i2goLiUlPqxtKwfw6GjpbhcEOf1UFJWzpHiMmKi3Oe9APL5FnUv4PfW2it9X/8GwFr7p9O9RkUtEjgLN+bx4qwNLN2y/8R9qQ1i6doikYuaxRMf7TlRPoXFpWTvO0L2vkKy9x/h+5yDlJRZ3C5D2UnXHomKcHGsrJyGcVFYC3sPF59TtsgIF8dKz77qjTEQ6DPpTx0jgMtAZZdccbsqWv/U55/pPeOiIig8Vnri/eKiImjbJI4Pf9YLc/z/IlVwpqL259OJZsD2k77OAS6p5JvcC9wL0KJFiyqHFBH/9GnVgD6tGrB21yEysg8wqG1DmtWN9uu1RSVlLNm8jyWb99MwLor0lERSG8QS5/VQXm5x+Qor/0gJa3cfYlNuAW5jiIxwUVxaTp2oCA4VlZB3+Bjl1mIMuIwhNiqCqAgXW/MKaZzgpUPTBIpKysg5eJTcw8W0aliHjk3jKbeWdbsLWLMrnyYJ0TRLjKa4pBy3y5B/tIQdB44SHx1BVISbhGgPTet6KSopZ+fBo+w8eBRjKvZkk+IqjjV3GYiJjGDb/iOUW0ubRnG0qB/D5txClm7ehzfSTefmdXG7DBv2FrA7v4j46AgSYyIpKi2n6FgZRSUVU0XWQkykmyYJ0aQ0iGXj3sPkHy2hdcM44qMjWL3zEHsOFVEvNoqikjJyDxcT740gPtrD0WNl5BYUc6y0/JxK+mz82aMeClxlrR3t+3o4cIm19sHTvUZ71CIiVXO+lzndASSf9HVz330iIlID/CnqZUBrY0yqMSYSuBX4e2BjiYjIcWedow5znH4AAAS/SURBVLbWlhpjHgT+RcXheROttasDnkxERAD/PkzEWvsP4B8BziIiIpXQ9ahFRIKcilpEJMipqEVEgpyKWkQkyAXk6nnGmFwg+xxf3gDIq8Y4TgmXcUD4jEXjCD7hMpbqGEdLa21SZQ8EpKjPhzEm43Rn54SScBkHhM9YNI7gEy5jCfQ4NPUhIhLkVNQiIkEuGIv6TacDVJNwGQeEz1g0juATLmMJ6DiCbo5aRER+KBj3qEVE5CQqahGRIBc0RW2MucoYs84Ys9EY82un81SVMWarMWalMSbLGJPhu6+eMWamMWaD7+9Ep3Oeyhgz0Riz1xiz6qT7Ks1tKrzo20bfG2O6Opf8x04zlt8bY3b4tkuWMeaakx77jW8s64wxVzqT+seMMcnGmNnGmDXGmNXGmLG++0Nqu5xhHKG4TbzGmG+NMSt8Y/mD7/5UY8xSX+ZpvktBY4yJ8n290fd4ynkFsNY6/oeKy6duAtKASGAF0N7pXFUcw1agwSn3PQ382nf718BfnM5ZSe7+QFdg1dlyA9cAXwIG6AksdTq/H2P5PfB4Jc9t7/t3FgWk+v79uZ0egy9bE6Cr73YcFYtLtw+17XKGcYTiNjFAHd9tD7DU97P+ELjVd//rwH2+2/cDr/tu3wpMO5/vHyx71D2AjdbazdbaY8AHwBCHM1WHIcBk3+3JwI0OZqmUtXYesP+Uu0+XewgwxVZYAtQ1xjSpmaRnd5qxnM4Q4ANrbbG1dguwkYp/h46z1u6y1mb6bh8G1lKxdmlIbZczjON0gnmbWGttge9Lj++PBQYBH/vuP3WbHN9WHwOXmfNYTDFYirqyBXTPtEGDkQW+MsYs9y30C9DIWrvLd3s30MiZaFV2utyhup0e9E0JTDxp+ikkxuL7lbkLFXtwIbtdThkHhOA2Mca4jTFZwF5gJhV7/AettaW+p5yc98RYfI/nA/XP9XsHS1GHg77W2q7A1cADxpj+Jz9oK34HCrljIUM190leAy4AOgO7gOecjeM/Y0wd4BPgYWvtoZMfC6XtUsk4QnKbWGvLrLWdqVg3tgfQtqa+d7AUdcgvoGut3eH7ey/wKRUbcs/xX0F9f+91LmGVnC53yG0na+0e339g5cBb/OdX6aAeizHGQ0W5TbXWTvfdHXLbpbJxhOo2Oc5aexCYDfSiYprp+EpZJ+c9MRbf4wnAvnP9nsFS1CG9gK4xJtYYE3f8NnAFsIqKMdzte9rdwGfOJKyy0+X+O3CX7yiDnkD+Sb+KB6VT5mpvomK7QMVYbvV9Op8KtAa+rel8lfHNZU4A1lprnz/poZDaLqcbR4hukyRjTF3f7WhgMBVz7rOBob6nnbpNjm+rocA3vt+Czo3Tn6ae9KnqNVR8KrwJ+K3TeaqYPY2KT6tXAKuP56diTmoWsAH4GqjndNZKsr9Pxa+fJVTMsY06XW4qPvl+xbeNVgLpTuf3Yyzv+LJ+7/uPp8lJz/+tbyzrgKudzn9Srr5UTGt8D2T5/lwTatvlDOMIxW3SCfjOl3kV8P9896dR8T+TjcBHQJTvfq/v642+x9PO5/vrFHIRkSAXLFMfIiJyGipqEZEgp6IWEQlyKmoRkSCnohYRCXIqahGRIKeiFhEJcv8fMo+y+IztoG4AAAAASUVORK5CYII=\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    }
  ]
}
	{
	"nbformat": 4,
	"nbformat_minor": 0,
	"metadata": {
	"colab": {
	"name": "Finding points for given distances.ipynb",
	"provenance": [],
	"collapsed_sections": [],
	"authorship_tag": "ABX9TyNyWO9Nsapp3Vz+azPgTTm8",
	"include_colab_link": true
	},
	"kernelspec": {
	"name": "python3",
	"display_name": "Python 3"
	},
	"language_info": {
	"name": "python"
	}
	},
	"cells": [
	{
	"cell_type": "markdown",
	"metadata": {
	"id": "view-in-github",
	"colab_type": "text"
	},
	"source": [
	"<a href=\"https://colab.research.google.com/gist/ishay2b/7aba5d8c7832fe63eedcffab33dd2714/finding-points-for-given-distances.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
	]
	},
	{
	"cell_type": "code",
	"execution_count": 2,
	"metadata": {
	"colab": {
	"base_uri": "https://localhost:8080/"
	},
	"id": "MR5udthOizcS",
	"outputId": "958fbd75-8168-4ad7-e706-a7b921843b63"
	},
	"outputs": [
	{
	"output_type": "stream",
	"name": "stdout",
	"text": [
	"Success reaching total distance 0.002 after 303 iterations\n",
	"Possible Found positions:\n",
	"[[3.9977682 0.52691495]\n",
	" [5.2908134 3.2339554 ]\n",
	" [0.38747454 2.249567 ]]\n",
	"Accumlated error 0.0013997554779052734\n",
	"Actual distnaces : [3.0000055 4.000219 5.0011754]\n",
	"Vs requested distnaces: [3. 4. 5.]\n"
	]
	}
	],
	"source": [
	"import torch\n",
	"from numpy import array as _A, dtype\n",
	"from torch.autograd import Variable\n",
	"import torch.optim as optim\n",
	"import numpy as np\n",
	"import matplotlib.pyplot as plt\n",
	"\n",
	"N = 3\n",
	"device = \"cpu\"\n",
	"\n",
	"\n",
	"def pairs(n):\n",
	" ''' Create all valid non repeating pairs '''\n",
	" return [[a, b] for a in range(n) for b in range(a, n) if a != b]\n",
	"\n",
	"\n",
	"all_pairs = _A(pairs(N))\n",
	"src_index = all_pairs[:, 0] # Index of the points p1\n",
	"trg_index = all_pairs[:, 1] # Index of the points p2\n",
	"\n",
	"#The input distances vector, assume ordered same as the pairs function above.\n",
	"all_distances = torch.tensor([3.0, 4.0, 5.0], requires_grad=False, device=device)\n",
	"\n",
	"#Random starting points scaled to mean distance as nicer starting points.\n",
	"start_points = np.random.random((3, 2)) * all_distances.mean().detach().cpu().numpy()\n",
	"\n",
	"# The coordintes we are trying to find.\n",
	"positions = torch.tensor(start_points, requires_grad=True,\n",
	" device=device, dtype=torch.float) \n",
	"\n",
	"optimizer = optim.Adam([positions], lr=0.01)\n",
	"\n",
	"loss_graph = [] # For visual\n",
	"threshold = 0.002 # Min accumulated error for early stopping.\n",
	"\n",
	"for i in range(1000):\n",
	" optimizer.zero_grad()\n",
	" p1 = positions[src_index] # Vector of start points.\n",
	" p2 = positions[trg_index] # Vector of end points.\n",
	" euc = ((p1 - p2)2).sum(axis=1)0.5 # Vector of euclidian distances.\n",
	" loss = (euc - all_distances).abs().sum() # Total loss error to minimize.\n",
	" loss_graph.append(float(loss.detach().cpu().numpy()))\n",
	" if loss_graph[-1] < threshold:\n",
	" print(f\"Success reaching total distance {threshold} after {i+1} iterations\")\n",
	" break\n",
	" loss.backward() # Now Optimize!\n",
	" optimizer.step()\n",
	"\n",
	"\n",
	"print(\"Possible Found positions:\")\n",
	"print(positions.detach().cpu().numpy())\n",
	"print(f\"Accumlated error {loss_graph[-1]}\")\n",
	"print(\"Actual distnaces :\", euc.detach().cpu().numpy())\n",
	"print(\"Vs requested distnaces:\", all_distances.detach().cpu().numpy())\n",
	"\n"
	]
	},
	{
	"cell_type": "code",
	"source": [
	"plt.plot(loss_graph, label=\"loss_graph\")\n",
	"plt.legend();"
	],
	"metadata": {
	"colab": {
	"base_uri": "https://localhost:8080/",
	"height": 265
	},
	"id": "1pm-HyUci7I9",
	"outputId": "b9031b33-6cac-4730-ebca-22e4a2a78ed1"
	},
	"execution_count": 3,
	"outputs": [
	{
	"output_type": "display_data",
	"data": {
	"image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAD4CAYAAADFAawfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3dd3hUZcL+8e8zk0kmCUkIEHogiSBVpASkgyh2xcIuNlSau7bFuuVd9/fuvuv77q51sRdAQFGx4LriuisivUmIQdrSCYSaUAIJJKQ8vz8ysIgBJpDJmZncn+viYjIt95ODtyfPnHMeY61FRESCl8vpACIicmYqahGRIKeiFhEJcipqEZEgp6IWEQlyEYF40wYNGtiUlJRAvLWISFhavnx5nrU2qbLHAlLUKSkpZGRkBOKtRUTCkjEm+3SPaepDRCTIqahFRIKcilpEJMj5NUdtjKkLjAc6AhYYaa1dHMhgIuKMkpIScnJyKCoqcjpKWPJ6vTRv3hyPx+P3a/z9MHEc8E9r7VBjTCQQcy4BRST45eTkEBcXR0pKCsYYp+OEFWst+/btIycnh9TUVL9fd9apD2NMAtAfmOD7RsestQfPOamIBLWioiLq16+vkg4AYwz169ev8m8r/sxRpwK5wNvGmO+MMeONMbGVBLjXGJNhjMnIzc2tUggRCS4q6cA5l5+tP0UdAXQFXrPWdgEKgV+f+iRr7ZvW2nRrbXpSUqXHbJ9RUUkZb87bxJLN+6r8WhGRcOZPUecAOdbapb6vP6aiuKs3iDFMWLCFl7/ZWN1vLSIS0s5a1Nba3cB2Y0wb312XAWuqO0hkhIu7e6ewYGMea3Yequ63F5EQUqdOHacjnLOUlBTy8vKq9T39PY76IWCqMeZ7oDPwf9WawueOHi2JiXQzfsHmQLy9iIjfSktLnY5wgl+H51lrs4D0AGchIcbDT9OTeXdJNr+8si2NE7yB/pYicgZ/+Hx1tf+G275pPP99fQe/nmut5Ze//CVffvklxhiefPJJhg0bxq5duxg2bBiHDh2itLSU1157jd69ezNq1CgyMjIwxjBy5EgeeeSRSt932bJljBo1CpfLxeDBg/nyyy9ZtWoVkyZNYvr06RQUFFBWVsYXX3zBkCFDOHDgACUlJTz11FMMGTKErVu3ctVVV9GtWzcyMzPp0KEDU6ZMISam4sjll156ic8//5ySkhI++ugj2rZte14/s6A7M3FU31TKrWXSoq1ORxERh02fPp2srCxWrFjB119/zRNPPMGuXbt47733uPLKK0881rlzZ7KystixYwerVq1i5cqVjBgx4rTvO2LECN544w2ysrJwu90/eCwzM5OPP/6YuXPn4vV6+fTTT8nMzGT27Nk89thjHF9ndt26ddx///2sXbuW+Ph4Xn311RPv0aBBAzIzM7nvvvt49tlnz/vnEJCr552P5HoxXNWxMe8tzeahQa2IjQq6iCK1hr97voGyYMECbrvtNtxuN40aNWLAgAEsW7aM7t27M3LkSEpKSrjxxhvp3LkzaWlpbN68mYceeohrr72WK664otL3PHjwIIcPH6ZXr14A3H777cyYMePE44MHD6ZevXpAxR79f/3XfzFv3jxcLhc7duxgz549ACQnJ9OnTx8A7rzzTl588UUef/xxAG6++WYAunXrxvTp08/75xB0e9QAo/ulcaiolA8ztjsdRUSCUP/+/Zk3bx7NmjXjnnvuYcqUKSQmJrJixQoGDhzI66+/zujRo8/pvWNj/3OayNSpU8nNzWX58uVkZWXRqFGjEyernHo89MlfR0VFAeB2u6tlrjsoi7pri0TSWyYyYcEWSsvKnY4jIg7p168f06ZNo6ysjNzcXObNm0ePHj3Izs6mUaNGjBkzhtGjR5OZmUleXh7l5eXccsstPPXUU2RmZlb6nnXr1iUuLo6lSyuOOP7ggw9O+/3z8/Np2LAhHo+H2bNnk539n0tGb9u2jcWLKy559N5779G3b99qHPkPBe28wuh+afz83eX8a/Ueru3UxOk4IuKAm266icWLF3PxxRdjjOHpp5+mcePGTJ48mWeeeQaPx0OdOnWYMmUKO3bsYMSIEZSXV+zc/elPfzrt+06YMIExY8bgcrkYMGAACQkJlT7vjjvu4Prrr+eiiy4iPT39Bx8KtmnThldeeYWRI0fSvn177rvvvuod/EnM8Ynx6pSenm7Pd4WXsnLLoOfmkBgTyaf399YprSI1ZO3atbRr187pGAFVUFBw4ljtP//5z+zatYtx48b5/fqtW7dy3XXXsWrVqnP6/pX9jI0xy621lR5dF5RTHwBul2FU31Syth9kefYBp+OISBj54osv6Ny5Mx07dmT+/Pk8+eSTTkc6o6Cd+gAY2q05z89cz1vzN5OeUs/pOCISYh544AEWLlz4g/vGjh3LiBEjGDZs2Dm/b0pKyjnvTZ+LoC7qmMgI7rykJa/M2ciWvEJSG/zoon0iEgDW2rCYbnzllVecjvAj5zLdHLRTH8fd1bslHpeLiQu2OB1FpFbwer3s27fvnApFzuz4wgFeb9XOug7qPWqAhnFehnRuykfLt/Po4AtJjI10OpJIWGvevDk5OTnouvKBcXwprqoI+qIGGNM/jY+W5zB1aTYPDmrtdByRsObxeKq0TJQEXtBPfQBc2CiOARcmMWlRNkUlZU7HERGpUSFR1ABj+qWRV1DM37N2Oh1FRKRGhUxR92lVn7aN4xi/YLM+5BCRWiVkitoYw73901i/p4C56/Uhh4jUHiFT1ADXdWpKo/go3pqvFWBEpPYIqaKOjHBxT+9UFm7cx+qd+U7HERGpESFV1AC392hBTKSbCfN1AoyI1A4hV9QJMR6GdU/m7yt2sju/yOk4IiIBF3JFDTCyj9ZVFJHaIySLOrleDFd3bMLUpdkUFAfPku4iIoEQkkUNMLpfKoeLSvlwmdZVFJHwFrJF3cW3ruLEhVpXUUTCm19FbYzZaoxZaYzJMsac3xpb1WhM/zRyDhzlX6v3OB1FRCRgqrJHfam1tvPp1vRywuXtGpFSP4Y35+u0chEJXyE79QH/WVdxxfaDZGhdRREJU/4WtQW+MsYsN8bcW9kTjDH3GmMyjDEZNXnB8aHdkqkb4+GteTqtXETCk79F3dda2xW4GnjAGNP/1CdYa9+01qZba9OTkpKqNeSZREe6Gd6zJTPX7mFLXmGNfV8RkZriV1Fba3f4/t4LfAr0CGSoqhreq2JdxQkLtFctIuHnrEVtjIk1xsQdvw1cAdTcOul+aBjn5cYuTfl4eQ4HCo85HUdEpFr5s0fdCFhgjFkBfAt8Ya39Z2BjVd3ofmkUlZTz7pJsp6OIiFSrsy5ua63dDFxcA1nOy4WN4hjYJonJi7MZ0z8Nr8ftdCQRkWoR0ofnnUrrKopIOAqrou59QX3aNYnnLZ0AIyJhJKyK2hjDmH6pbNhbwBytqygiYSKsihoq1lVsHO9lvNZVFJEwEXZFHRnh4p4+KVpXUUTCRtgVNcBtPVoQG+lmvNZVFJEwEJZFnRDt4afdk/l8xU525R91Oo6IyHkJy6IGrasoIuEjbIs6uV4MV1/UhPeWbtO6iiIS0sK2qKHiBJjDRaVM07qKIhLCwrqoOyfXpXtKIhMXaF1FEQldYV3UULFXvePgUf65erfTUUREzknYF/Xl7RqR2iCWt+Zv0WnlIhKSwr6oXS7DSK2rKCIhLOyLGmBo1+Ykxnh4U+sqikgIqhVFHR3p5s6eLfla6yqKSAiqFUUNcFevFK2rKCIhqdYUdVJcFDd1acbHy3PYr3UVRSSE1JqiBhjdL1XrKopIyKlVRd3at67ilMVbKSopczqOiIhfalVRw/F1FY9pXUURCRm1rqh7X1Cfto3jGL9A6yqKSGiodUVdsa5iGuv3FDB/Q57TcUREzqrWFTXA9Rc3pWFcFG9pXUURCQG1sqgjI1zc3TuF+RvyWLf7sNNxRETOyO+iNsa4jTHfGWNmBDJQTbnjkhZEe9xarVxEgl5V9qjHAmsDFaSm1Y2JZGi35nyWtZO9h4ucjiMiclp+FbUxpjlwLTA+sHFq1si+qZSUl/PuYp0AIyLBy9896r8CvwROu0yKMeZeY0yGMSYjNze3WsIFWmqDWC5v14h3lmRz9JhOgBGR4HTWojbGXAfstdYuP9PzrLVvWmvTrbXpSUlJ1RYw0Mb0S+PAkRKmf5fjdBQRkUr5s0fdB7jBGLMV+AAYZIx5N6CpalD3lEQ6NU9gwvwtlJfrBBgRCT5nLWpr7W+stc2ttSnArcA31to7A56shhhjGNU3lc15hcxet9fpOCIiP1Irj6M+1TUXNaFpgpfx87c4HUVE5EeqVNTW2jnW2usCFcYpHreLe/qksHjzPlbtyHc6jojID2iP2mdY9xbERrqZsEB71SISXFTUPgnRHoZ1b8HnK3ayK/+o03FERE5QUZ9kRJ8Uyq1l8iKdACMiwUNFfZLkejFc3bEJ7y3NprC41Ok4IiKAivpHRvVL5VBRKR9lbHc6iogIoKL+ka4tEunaoi4TF26lTCfAiEgQUFFXYky/NLbtP8LMNbudjiIioqKuzBUdGpNcL1onwIhIUFBRV8LtMozsk0pG9gG+23bA6TgiUsupqE/jJ+nJxHkjGK8TYETEYSrq06gTFcHtPVrw5cpdbN9/xOk4IlKLqajP4J4+KbiMYdKirU5HEZFaTEV9Bk0Sorm2UxOmLdvOoaISp+OISC2loj6L0X3TKCgu5cNlOgFGRJyhoj6Li5oncElqPd5euJXSstMuGSkiEjAqaj+M7pfGjoNH+XKVToARkZqnovbDZW0bktoglvHzN2OtTisXkZqlovaDy2UY2TeVFTn5ZGTrBBgRqVkqaj8N7dqcujEexs/f7HQUEallVNR+io50c+clLflqzR625hU6HUdEahEVdRXc1aslHpeLtxfqtHIRqTkq6ipoGO/lhs5N+TAjh4NHjjkdR0RqCRV1FY3qm8rRkjLe+3ab01FEpJZQUVdRuybx9G3VgMmLtnKsVCfAiEjgnbWojTFeY8y3xpgVxpjVxpg/1ESwYDa6Xyp7DhUz4/udTkcRkVrAnz3qYmCQtfZioDNwlTGmZ2BjBbcBFybRumEdxs/fohNgRCTgzlrUtkKB70uP70+tbidjDKP7pbJm1yEWb9rndBwRCXN+zVEbY9zGmCxgLzDTWru0kufca4zJMMZk5ObmVnfOoDOkczPqx0ZqBRgRCTi/itpaW2at7Qw0B3oYYzpW8pw3rbXp1tr0pKSk6s4ZdLweN8N7teSbf+9l496Cs79AROQcVemoD2vtQWA2cFVg4oSW4T1bEhnh4q15Oq1cRALHn6M+kowxdX23o4HBwL8DHSwU1K8TxW3dk5mWsZ1HpmXp1HIRCYgIP57TBJhsjHFTUewfWmtnBDZW6Pjdde2Jj/bwxrzNfJa1g+s6NeWJK9uQXC/G6WgiEiZMIA4vS09PtxkZGdX+vsFs7+EiJi7YyqRFWyi38PP+aTwwqBVREW6no4lICDDGLLfWplf2mM5MrCYN47z8+uq2fPPYQK7u2JgXv9nI9S8tYGVOvtPRRCTEqairWdO60Yy7tQtvj+hO/tESbnx1Ic/PXK/1FkXknKmoA+TSNg356uEB3HBxU16ctYG7Jn5LXkGx07FEJASpqAMoIcbDC8M68/QtncjIPsD1Ly0gc5uW8hKRqlFR14Cfdk9m+n29iXAbhr2xmPHzN7Pz4FGnY4lIiNBRHzUo/0gJD76fyfwNeQD0vzCJX17Zho7NEhxOJiJO01EfQSIhxsOUkT2Y8VBfHh18Iat25DPklYW8NW+zrsInIqeloq5hxhg6NkvgF5e1ZvbjAxncrhH/+4+13D81k6KSMqfjiUgQUlE7KCHaw2t3duW317Tjn6t3M2ryMo4cK3U6logEGRW1w4wxjOmfxvM/vZjFm/Zx98RvOVxU4nQsEQkiKuogcVOX5rx0W1e+23aQO8YvZfv+I05HEpEgoaIOItd2asLrd3Zjc24hV4+bz8fLc5yOJCJBQEUdZC5v34gvx/ajfdN4Hv9oBb+ZvlKrnYvUcirqIJRcL4b3x/TkvoEX8P6327j9rSXsPVzkdCwRcYiKOki5XYZfXdWWl2/vwuqdh7jhpYWs233Y6Vgi4gAVdZC7rlNTPrmvN+XW8tM3FrM8W9cKEaltVNQhoH3TeD65rzeJMR7uHL+UuevDf5V3EfkPFXWISK4Xw0c/701qg1hGT17GV6t3Ox1JRGqIijqEJMVF8cHPetKhaQL3T81UWYvUEirqEBPv9TBlVA86NKso60kLt1Bergs6iYQzFXUIivd6eGdUDwZcmMTvP1/DtS8t4Lmv1rFhj44KEQlHKuoQFe/18NZd6Tz3k4uJdBtenbOJwS/M4xfvf8cOLUogEla0cECY2FdQzKRFW3lz3mYA/mdIB4Z1b+FwKhHxlxYOqAXq14nisSva8M3jA0lPSeRXn6zkiY9WcPSYrnEtEurOWtTGmGRjzGxjzBpjzGpjzNiaCCbnplndaKaMvIRfDGrFx5k53PTqQrL3FTodS0TOgz971KXAY9ba9kBP4AFjTPvAxpLz4XYZHr2iDW/f051d+UXc8PJCFvjWaRSR0HPWorbW7rLWZvpuHwbWAs0CHUzO38A2Dfn7g31oHO9l+MSljP3gOzK27tf6jCIhpkofJhpjUoB5QEdr7aFTHrsXuBegRYsW3bKzs6svpZyXwuJSXpy1gXeXZFN4rIxuLRP5/fUduKi5Vj8XCRZn+jDR76I2xtQB5gL/a62dfqbn6qiP4FRYXMr073bw0qwN7Cs8xph+aTx8eWu8HrfT0URqvfM+6sMY4wE+AaaeraQleMVGRTC8Z0tmPjqAoV2b8/rcTVw9bj5LN+9zOpqInIE/R30YYAKw1lr7fOAjSaAlRHv4y9BOTB19CaXl5Qx7cwlP/m0lhcVaAV0kGPmzR90HGA4MMsZk+f5cE+BcUgP6tGrAvx7uz+i+qUxduo3rXlrAqh35TscSkVPozEQBYOnmfTw8LYu8gmKeuLINI/ukEuHW+VAiNUVnJspZXZJWny/H9uOyto34v3/8m75/mc24rzeQf7TE6WgitZ6KWk6oGxPJa3d2Zfxd6VzYOI4Xvl7PoGfnMG3ZNl1KVcRBKmr5AWMMl7dvxJSRPZjxUF9SG8Tyq09WcuOrC1m9U/PXIk5QUctpdWyWwEc/78Vfh3VmV34RQ15eyLivN1BSVu50NJFaRUUtZ2SM4cYuzfjq4f5c26kJL3y9njvGL9U1r0VqkIpa/JIYG8m4W7vwwrCLWZmTz2XPzeGFmespKtFlVEUCTUUtVXJTl+Z8/dgALm/XiHGzNjD4hbnMXrfX6VgiYU1FLVXWrG40L9/elffH9CTS7WLE28u4f+pyducXOR1NJCypqOWc9bqgPl+O7c8TV7Zh1tq9XPbcHCYu2KIPG0WqmYpazktkhIsHLm3FzEcGkJ5Sj/+ZsYZef5rFhAVbKFVhi1QLFbVUixb1Y5g0ojuTR/agXZN4/jhjDcMnfKvpEJFqoKKWamOMYcCFSUwZ2YNnhnZi+bYD9H9mNn/4fDV7D6uwRc6VilqqnTGGn6QnM+vRAdzUuRlTFmfT/+mKa4ccK9V0iEhVqaglYJLrxfCXoZ34xnc43wtfr+f6lxbw3bYDTkcTCSkqagm4lvVjefn2rky4O538oyXc/Noi/jhjDUeOaaECEX+oqKXGXNauETMf7c8dl7RgwoItDH5+Hv9ctVuroouchYpaalSc18NTN17Exz/vRZw3gp+/u5wRk5aRva/Q6WgiQUtFLY5IT6nHjIf68uS17Vi2ZT+DX5jH81+t07qNIpVQUYtjItwuRvdLY9ZjA7mqQ2Ne/GYjlz47h799t0PTISInUVGL4xoneHnxti58cl9vmiR4eXhaFndN/JZNuQVORxMJCipqCRrdWiYy/f4+/HFIB7K2HeSKF+bxm+nfs+eQTpaR2k1FLUHF7TIM75XC7CcGMrxnSz5ensPlz83lnSXZWrdRai0VtQSlBnWi+P0NHZj5yAA6JSfwu7+t4pbXF2ndRqmVVNQS1FIaxPLuqEt4/qcXk73vCNe9tIBHp2WxU0uBSS1y1qI2xkw0xuw1xqyqiUAipzLGcHPX5sx+bCD39k9jxspdDH5+LuPnb6a4VEuBSfjzZ496EnBVgHOInFVCjIffXN2OWY9WXPv6qS/Wctlzc/l8xU4dzidh7axFba2dB+yvgSwifkmuF8PkkT14Z1QP4rweHnr/O37y+mKWZ+tiTxKeNEctIatf6yRmPNSXv9xyEVv3FXLLa4u4Y/wSMnV1Pgkz1VbUxph7jTEZxpiM3Nzc6npbkTNyuwzDurdg7hOX8ttr2rFudwG3vLaI3366kvyjJU7HE6kWxp+5PWNMCjDDWtvRnzdNT0+3GRkZ55dM5BwUFpfy/Mz1vL1wC/XrRDH2stYM6dyUOK/H6WgiZ2SMWW6tTa/sMU19SFiJjYrgd9e157MH+tIkwcuTf1vFpc/O4bMsXT9EQpc/h+e9DywG2hhjcowxowIfS+T8XNQ8gc8e6MMn9/Wmad1oxn6QxfAJ37J+z2Gno4lUmV9TH1WlqQ8JJmXllncWb+X5mespPFbGHZe04BeXtaZBnSino4mcoKkPqdXcLsM9fVKZ88Sl3NYjmalLtzHwmTm8/M0Gjh7TCTMS/FTUUmvUi43kqRsv4l8P96fXBfV59qv1mr+WkKCillqnVcM6vHVXOh/+rBcN46MY+0EWoydnsCtf1w+R4KSillqrR2o9Pr2/D7+7rj0LN+Ux4Jk5/PdnqzhQeMzpaCI/oKKWWs3tMozqm8rMRwZwS9dmvLt0GwOfncOkhVsoKSt3Op4IoKM+RH5g3e7D/HHGGhZszCO1QSwPXNqKIZ2b4nFrn0YCS0d9iPipTeM43hnVg/F3pRPtcfP4Ryu48q/zWLJ5n9PRpBZTUYucwhjD5e0b8cUv+vLm8G6UlJVz65tLeHRaFlvyCp2OJ7VQhNMBRIKVMYYrOjSmX+skxs3awNsLt/D3FTu5q1cKYy9rTUKMrh8iNUN71CJnER3p5tdXt2XBrwbxk/Rk3l60hYHPzuadJdmU6gNHqQEqahE/JcVF8aebL+KLh/rRpnEcv/vbKq4eN5/Z6/Y6HU3CnIpapIraN43n/TE9ecM3fz3i7WXcPfFbtmr+WgJERS1yDowxXNmhMV89MoAnr23H8uwDXP78XH776Uqd4SjVTsdRi1SDvYeKePGbDUxbth2D4Y6eLfjFoNYkxkY6HU1CxJmOo1ZRi1Sj7fuP8PI3G/lo+XbivB5+cVlrhvdsSWSEfnmVM9MJLyI1JLleDH8Z2okvx/anU/ME/jhjDVf+dR4zvt9Jebmu0CfnRkUtEgBtGscxZWQP3r6nO5FuFw++9x03vLKAuetzdUlVqTIVtUiAGGO4tG1D/jG2Hy8Mu5iDR0q4e+K33PbWEjK3HXA6noQQFbVIgLldhpu6NOebxwbyhxs6sHFvATe/uogxUzK0hqP4RR8mitSwwuJS3l64hTfmbqbgWCk3d2nOw5e3JrlejNPRxEE66kMkCB0oPMZrczcxadFWrLXccUlLHhzUSovu1lIqapEgtiv/KC/O2sCHGTlERbgY3S+NMf1SifPqok+1iYpaJARsyi3g+a/W88XKXSTGeLh/YCvu7NmS6Ei309GkBqioRULIypx8nv7Xv5m/IY8GdSIZ0y+NO3u2JDZKVyUOZypqkRC0bOt+Xpy1gfkb8kiM8TC6Xxp39WqpKZEwdd5FbYy5ChgHuIHx1to/n+n5KmqR6pO57QAvzdrA7HW51I3x8MDAVgzv1RKv54dTItZacguKOVJcRlFpGSWlluhIN80To3/0XAk+51XUxhg3sB4YDOQAy4DbrLVrTvcaFbVI9fs+5yDPfrWeeetzaZrg5e7eKfS/MIl563PJyD5AZvYB9hUe+9HrjIEm8V6aJUbTtUUiPdPq07tVfQ4UlrA8+wA5B46QV1BMVISb3YeKMEBaUh3aNYnjwkZxGAPFJeVERrg4XFRKfHQEDepEEeEyWAvl1lJYXMaRklIaxnnZX3iMtbsOERvlplndGOrFRrJ21yE27i3A7TJckFSHCxvXYU9+MbkFxXg9LkrLLAnRHprU9XLoaCml5eUcPFLCzoNH8XrcNK0bTZMELy5jOFxUwp5DxRhTMb6C4lIaxXmJ8rhYtSOfPYeKaZLgpXtqPQqLS9mUW4C1kJwYQ3K9aA4eKaGguBSvx43X48LrceNxu7DWUlJm2ZRbwO5DRTSrG03DuCg27i2goLiUlPqxtKwfw6GjpbhcEOf1UFJWzpHiMmKi3Oe9APL5FnUv4PfW2it9X/8GwFr7p9O9RkUtEjgLN+bx4qwNLN2y/8R9qQ1i6doikYuaxRMf7TlRPoXFpWTvO0L2vkKy9x/h+5yDlJRZ3C5D2UnXHomKcHGsrJyGcVFYC3sPF59TtsgIF8dKz77qjTEQ6DPpTx0jgMtAZZdccbsqWv/U55/pPeOiIig8Vnri/eKiImjbJI4Pf9YLc/z/IlVwpqL259OJZsD2k77OAS6p5JvcC9wL0KJFiyqHFBH/9GnVgD6tGrB21yEysg8wqG1DmtWN9uu1RSVlLNm8jyWb99MwLor0lERSG8QS5/VQXm5x+Qor/0gJa3cfYlNuAW5jiIxwUVxaTp2oCA4VlZB3+Bjl1mIMuIwhNiqCqAgXW/MKaZzgpUPTBIpKysg5eJTcw8W0aliHjk3jKbeWdbsLWLMrnyYJ0TRLjKa4pBy3y5B/tIQdB44SHx1BVISbhGgPTet6KSopZ+fBo+w8eBRjKvZkk+IqjjV3GYiJjGDb/iOUW0ubRnG0qB/D5txClm7ehzfSTefmdXG7DBv2FrA7v4j46AgSYyIpKi2n6FgZRSUVU0XWQkykmyYJ0aQ0iGXj3sPkHy2hdcM44qMjWL3zEHsOFVEvNoqikjJyDxcT740gPtrD0WNl5BYUc6y0/JxK+mz82aMeClxlrR3t+3o4cIm19sHTvUZ71CIiVXO+lzndASSf9HVz330iIlID/CnqZUBrY0yqMSYSuBX4e2BjiYjIcWedow5znH4AAAS/SURBVLbWlhpjHgT+RcXheROttasDnkxERAD/PkzEWvsP4B8BziIiIpXQ9ahFRIKcilpEJMipqEVEgpyKWkQkyAXk6nnGmFwg+xxf3gDIq8Y4TgmXcUD4jEXjCD7hMpbqGEdLa21SZQ8EpKjPhzEm43Rn54SScBkHhM9YNI7gEy5jCfQ4NPUhIhLkVNQiIkEuGIv6TacDVJNwGQeEz1g0juATLmMJ6DiCbo5aRER+KBj3qEVE5CQqahGRIBc0RW2MucoYs84Ys9EY82un81SVMWarMWalMSbLGJPhu6+eMWamMWaD7+9Ep3Oeyhgz0Riz1xiz6qT7Ks1tKrzo20bfG2O6Opf8x04zlt8bY3b4tkuWMeaakx77jW8s64wxVzqT+seMMcnGmNnGmDXGmNXGmLG++0Nqu5xhHKG4TbzGmG+NMSt8Y/mD7/5UY8xSX+ZpvktBY4yJ8n290fd4ynkFsNY6/oeKy6duAtKASGAF0N7pXFUcw1agwSn3PQ382nf718BfnM5ZSe7+QFdg1dlyA9cAXwIG6AksdTq/H2P5PfB4Jc9t7/t3FgWk+v79uZ0egy9bE6Cr73YcFYtLtw+17XKGcYTiNjFAHd9tD7DU97P+ELjVd//rwH2+2/cDr/tu3wpMO5/vHyx71D2AjdbazdbaY8AHwBCHM1WHIcBk3+3JwI0OZqmUtXYesP+Uu0+XewgwxVZYAtQ1xjSpmaRnd5qxnM4Q4ANrbbG1dguwkYp/h46z1u6y1mb6bh8G1lKxdmlIbZczjON0gnmbWGttge9Lj++PBQYBH/vuP3WbHN9WHwOXmfNYTDFYirqyBXTPtEGDkQW+MsYs9y30C9DIWrvLd3s30MiZaFV2utyhup0e9E0JTDxp+ikkxuL7lbkLFXtwIbtdThkHhOA2Mca4jTFZwF5gJhV7/AettaW+p5yc98RYfI/nA/XP9XsHS1GHg77W2q7A1cADxpj+Jz9oK34HCrljIUM190leAy4AOgO7gOecjeM/Y0wd4BPgYWvtoZMfC6XtUsk4QnKbWGvLrLWdqVg3tgfQtqa+d7AUdcgvoGut3eH7ey/wKRUbcs/xX0F9f+91LmGVnC53yG0na+0e339g5cBb/OdX6aAeizHGQ0W5TbXWTvfdHXLbpbJxhOo2Oc5aexCYDfSiYprp+EpZJ+c9MRbf4wnAvnP9nsFS1CG9gK4xJtYYE3f8NnAFsIqKMdzte9rdwGfOJKyy0+X+O3CX7yiDnkD+Sb+KB6VT5mpvomK7QMVYbvV9Op8KtAa+rel8lfHNZU4A1lprnz/poZDaLqcbR4hukyRjTF3f7WhgMBVz7rOBob6nnbpNjm+rocA3vt+Czo3Tn6ae9KnqNVR8KrwJ+K3TeaqYPY2KT6tXAKuP56diTmoWsAH4GqjndNZKsr9Pxa+fJVTMsY06XW4qPvl+xbeNVgLpTuf3Yyzv+LJ+7/uPp8lJz/+tbyzrgKudzn9Srr5UTGt8D2T5/lwTatvlDOMIxW3SCfjOl3kV8P9896dR8T+TjcBHQJTvfq/v642+x9PO5/vrFHIRkSAXLFMfIiJyGipqEZEgp6IWEQlyKmoRkSCnohYRCXIqahGRIKeiFhEJcv8fMo+y+IztoG4AAAAASUVORK5CYII=\n",
	"text/plain": [
	"<Figure size 432x288 with 1 Axes>"
	]
	},
	"metadata": {
	"needs_background": "light"
	}
	}
	]
	}
	]
	}