{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "This code is meant to simulate the orbital evolution of Leela's constellation from Futurama episode \"Time Keeps on Slippin'\". \n",
    "\n",
    "Keep in mind the orbit simulator <a href=\"http://rebound.readthedocs.io/en/latest/\">REBOUND</a> is really designed to evolve well-defined Keplerian orbits, which the stars in Leela's constellation must definitely do not occupy."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "import rebound\n",
    "import numpy as np\n",
    "import ipywidgets"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "First, grab image of Leela's constellation. I then made a png in gimp with one pixel defining each stellar center. That's what displayed below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgIAAAFeCAYAAAD34XfdAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAIABJREFUeJzt3Xl8VNXdP/DPmSUzmSVkDyQMSUCQAEa2EBZBRQSrRRYt\nvtCi/nysu9ZatWjVwqPVPq19Wh+wWrugL6lr3Z6iVmUThCAg8AQCGCAhCyQh+z7Z5vz+yKAhZCYz\nyczcmbmf9+t1XmHuPXPvNyeXzDfnnnuOkFKCiIiI1EmjdABERESkHCYCREREKsZEgIiISMWYCBAR\nEakYEwEiIiIVYyJARESkYkwEiIiIVIyJABERkYoxESAiIlIxRRIBIcQSIcRuIcQ2IcRWIcQ4JeIg\nIiJSO12gTyiEmAbgVQCTpZQnhBArAHwmhMiQUjYFOh4iIiI1U6JHYCWADVLKE87X69GdkNyqQCxE\nRESqFvAeAQBzAaw++0JKKYUQ3wCYB2CtEEILYHSv99QA4OpIRERE7gkAsb22HZNSdrl6Q0ATASFE\nHIAoABW9dlUAmOr892gARwIZFxERURjLAHDU1c5A3xowOb+29dre1mMfERERBUigE4EW51dDr+0G\nAM0BjoWIiEj1ApoISCmrAdQDSOq1ayiAgkDGQkRERMo8NbAZ348HgBBCAJgMYKNzU40CMREREYUr\nt5+rSiQCvwFwjRBilPP1TQA6ALzmfM2nA4iIiHzH7edqwB8flFLuEULcCuAtIUQrgC4AC6SUHCNA\nREQUYELK4PoDXAiRAOCM0nEQERGFiUQpZaWrnVx0iIiISMWYCBAREakYEwEiIiIVYyJARESkYkwE\niIiIVIyJABERkYoxESAiIlIxJgJEREQqxkSAiIhIxZgIEBERqRgTASIiIhUL+KJDROQ/VqsVS5cu\nhUajwfvvv4/6+nqvj5Geno4lS5bg2LFj+Ne//uWHKIkomLBHgCiMREZGYtasWbjkkktgMpkGdIz4\n+HgsWLAAEydOdFvPYDAgJSUF8fHxAzoPEQUHrj5IFEYMBgMyMjIghMCRI0dgt9u9PkZ0dDQyMjJQ\nU1ODb7/91mW9MWPG4Kc//SmOHTuGF154AcH2u4SIvuN29UHeGiAKI21tbThw4MCgjlFXV4ecnJx+\n63V1daG5uXlAyQYRBQ/2CBDRgOh0OlgsFnR2dqKpqUnpcIjINfYIEJHvdXZ2oq6uTukwiGiQOFiQ\niPpktVphs9lgNpuVDoWI/IiJABH1ae7cuVizZg1mzJihdChE5EdMBIioT11dXejo6IDD4VA6FCLy\nIw4WJKI+mUwmWK1WNDQ0oLW1VelwiGjgOFiQiLzX0tKClpYWpcMgIj/jrQEiIiIVYyJARESkYkwE\niIiIVIyJABERkYoxESCiPk2bNg0rV65EZmam0qEQkR8xESCiPl188cW47777MG7cOKVDISI/4uOD\nRNSnzZs348yZM4NezZCIghsnFCIiIgpvbicU4q0BIiIiFWMiQEREpGJMBIiIiFSMiQAREZGK8akB\nInJJo9Fg7NixMJlMOHr0KJqampQOiYh8jD0CROSSTqfDzTffjEcffRQpKSlKh0NEfsAeASJyyeFw\nYP/+/SgtLUVDQ4PS4RCRH3AeASIiovDGeQSIiIiob0wEiIiIVIyJABERkYoxESAiIlIxPjVARG5F\nRERACIGOjg44HA6lwyEiH2OPABG5xHkEiMIfEwEickkIgaSkJAwfPhwRERFKh0NEfsB5BIjIJSEE\nUlJSYDAYcOrUKdjtdqVDIiLvuZ1HgGMEiMglKSVKS0uVDoOI/Ii3BoiIiFSMiQAREZGK8dYAEbkk\nhEB6ejqMRiNOnjyJlpYWpUMiIh9jjwARuaTX63HLLbfw8UGiMMYeAaIwZzQaMXv2bERERGDbtm1o\nbGz0+L0OhwMFBQVobm5mbwBRmGIiQBTmrFYr7rzzTlitVhw+fNirRKCzsxOvvfaaV+cTQgDofuKA\niIIfEwGiMNfY2Ig///nPiIiIQFVVlV/PZbPZcNNNN+H06dN444030NnZ6dfzEdHgMREgCnN2ux1f\nfPFFQM4VGRmJjIwMGAyG73oGiCi4cWZBIvIZk8mE1NRU2O12nDx5krcHiIKD25kFmQgQERGFN7eJ\nAB8fJCIiUjEmAkRERCrGRICIiEjFmAgQERGpGBMBIiIiFWMiQEREpGJMBIiIiFSMiQAREZGKMREg\nIiJSMSYCREREKsZEgIiISMWYCBAREakYEwEiIiIVYyJARESkYkwEiIiIVEznqwMJIVYBWASgrsfm\nGinldT3q3AngJwDsznp3SClP+yoGIiIi8o7PEgEAEsBPpZTb+tophFgK4CkAF0kpa4QQTwLYIISY\nIqWUPoyDiIiIPOTrWwPCzb4nALwqpaxxvn4BwAQA1/g4BiIiIvJQQMYICCFiAUwEsPfsNillA4B8\nAPMCEQMRERGdz9eJwG1CiC1CiK+EEK8KIUY6t6c7v1b0ql/eYx8REREFmC8TgWIA+wBcAWA2gEIA\n3wghkgGYnHXaer2nHUCkD2MgIiIiL/gsEZBSrpNSviCldDgH/z2N7qcD7gHQ7Kxm6PU2A4AWX8VA\nRERE3vHbGAEppQPASQAjARQ4Nyf1qja0xz4iIiIKMJ8lAkKIF/rYnAKgWEpZB2A/gKk96kcBGA1g\no69iICIiIu/4skdgoRBi4dkXQojbAcQB+Ltz0zMAbnE+QQAADwA4KKX8xIcxEBERkRd8OaHQLwE8\nKIR4CEAEuscHzJNS5gOAlPIDIUQigC+EEHYANQAWujwaERER+Z0Itkn9hBAJAM4oHQcREVGYSJRS\nVrrayUWHiIiIVIyJABERkYoxESAiIlIxJgJEREQqxkSAiIhIxZgIEBH5mNlsRkxMDPR6vdKhEPWL\niQARkY8tW7YMq1atwujRo5UOhahfTASIiHzMYDDAbDZDp/PlnG1E/sEJhYiIfCw2NhYmkwlVVVWw\n2+1Kh0PkdkIhJgJE5BPz589HcnIyPvvsM5SVlXn9/ujoaCxYsABtbW347LPP0Nra6ocoiVTJbSLA\nfisi8okFCxYgKysLBw8ePCcR0Gq1iIyMhMPhQGtrK1z98RETE4Ply5ejvr4e27ZtYyLgASEETCYT\nhBBoaWmBw+FQOiQKQewRICKfmDJlCuLi4vDNN9+gurr6u+02mw233347zpw5g7/+9a9oa2vr8/0W\niwVTp05Fe3s79u7di/b29kCFHrIsFgtuv/12WCwW/OUvf0FFRYXSIVFwYo8AEfnfN9980+d2nU6H\n6Ojofu+VNzU1YevWrX6ILHhERkZi7Nix6OjowLfffouOjo5BHU8IgaioKFitVmi1Wh9FSWrDHgEi\n8iuDwYDExER0dHSgoqLC5a0BNUhNTcXq1atRX1+Pp556CvX19YM6nlarRWJiIrRaLSoqKgadWFDY\nYo8AEXlv8uTJmD17NjZu3Ii8vLwBH6etrQ0lJSU+jCx0NTc3Y9euXWhpafHJh3ZXV9eABmaGiuzs\nbEyZMgUbN25Efn6+T44phMC8efOQlJSEL774grdTwESAiFzIzs7Gk08+iaqqqkElAvS9qqoqvPzy\ny4M6hlarhdFoREdHh9/GUWg0GhgMBjgcjvPGdOh0OhiNRrS1tfm9B2LOnDm49957cfr0aZ8mAtde\ney0mTpyIvLw8JgJgIkBELmzevBm1tbXYvXu30qFQD+PHj8dtt92GnTt34p133vHLOUaMGIHbbrsN\nxcXFWLduHbq6ur7bN2PGDNx2221466238Nlnn/nl/Gd9+umnKCwsxL59+3x2TCkl1q9fj08++QRF\nRUU+O24oYyJARH369ttv8e233yodRlhKSEjABRdcgNLSUq9vmxiNRiQnJ2PIkCGDikEIgdGjR8Nk\nMuH48eNoamr6bl9ERASSkpLO2XaW2WzG8OHDYbVaB3V+Txw6dAiHDh3y6TGllPj66699esyQJ6UM\nqgIgAYBkYWFhCdeyYMECuXPnTvmTn/zE6/daLBY5evRomZiYOKgYdDqdfOaZZ+Sbb74pMzIyztln\nNBplenq6HDZsmBRCnLMvOjpaZmRkyNjYWMXbkcXjkuDuc5c9AkREAVZRUYGtW7eisLDQ6/c2NTXh\n2LFjg45BSom8vDxUV1ejsbHxnH12u91lbHV1dairqxv0+T2RlZWFSZMmYfPmzTh+/HhAzqlGTASI\niALswIEDOHDggMv9QojvljDu6OjwyyOXXV1dePPNN31+XF+aO3cu7r//flRWVjIR8CMmAkREQSY6\nOhq33norpJR47bXXUFtbq3RIivjkk09w8uRJnw4WpPMxESCiPg0bNgwjRoxAQUEBKitdzkVCfqDT\n6ZCcnPzdv9Xq4MGDOHjwoF+ObTabMXLkSDQ3N6OgoMAv5wgVGqUDIKLgNHfuXKxduxYzZsxQOhTV\nqaurw4svvogXX3wxYPfj1SY1NRWPP/44VqxYAY1G3R+F6k01icitU6dOYceOHThzhjN+B1pHRwdO\nnjypdBhemzlzJtLT07F161acOnVK6XDcampqwv79+1FaWqrqaa8BrjVA5JIQAhqNBlLKoFveVaPR\nQKPRoKurS/W/xGjgNBoNhBBwOBw+uY7+67/+C9deey3uuecebNmyxQcRko9wrQGigRgzZgxWrFiB\n3Nxcv83gNlCLFi3CjBkz8Prrr/vtHiqFN61WixUrViA5ORnr169HcXHxoI/53nvvYe/evZyIKsQw\nESBywWKx4MILL0RVVZXSoZxn6NChyMjICMjsbhR4ERERSE1NBQAUFRX5ZU0BIQRsNhtGjhyJyMhI\nnxxz9+7dfpuSOioqCmlpaaipqUFpaalfzqFWvDVA5ILVaoXNZkNDQ0PQ/eJJTk5GbGwsiouL0dDQ\noHQ45GOJiYl4/PHHIaXEc88955dxGkIIjBgxAkajESUlJWhpafH5OXwpOzsbv/zlL7Fp0ya88MIL\nSocTanhrgGggGhsbcfjwYaXD6NPp06dx+vRppcM4x+zZs5GYmIjt27dzgOEgtbe3fzd7oL9WGJRS\nhtSiOw0NDcjNzeWS1n7AHgEi8ok//elPmD59Ou6++24u6uIDQggA4GDQHtgmA8YeASLyvzfffBOb\nN28e0Pz5dL5Q/LC77rrrcPHFF+Mf//iHXwYM+qJNhBC46aabYLPZ8MYbb4RUr4i/MBGggIiPj4fN\nZkNxcTGqq6v7rX92mdOSkpI+l0Il/0tISEBycjJKSkpQU1PTb/3t27d7dNyUlBRERUXxZxuGRowY\ngczMTPzrX/9SOhSXhBBIT0/H2LFjYTablQ4nOCi97DCXIVZHWbp0qdy1a5dcvHixR/Ufeugh+dZb\nb8lJkyYpHrtayw033CC/+uorec011/j0uA8++KB899135eTJkxX/Hll8W2w2m8zMzJRRUVGKx+Ku\npKeny/Hjx0uz2ax4LAEqXIaYlFdTU4O8vDyPF08pKyvD8ePH0dzc7Lbe1KlTkZ6ejp07dwZ0JrNx\n48YhKysLO3fu9MmSsJ5KS0vDzJkzkZeXh//7v//z67mqq6uRl5fn8yluy8vLkZ+f3+/PlgIjISEB\nM2fORGVlJXbu3DmoY5WUlITEYD7evupF6R4A9gioo2g0GqnX66VGo/GovlarlTqdTgoh3Nb79a9/\nLfPz8+WVV14Z0O/n3nvvlU1NTQH/q3bJkiXy9OnTcuXKlUH3M/O0ePqzZQlMmTZtmty9e7dcu3at\n4rGw+K247RHQrlq1CsFk9erVZgCPKB0H+dbZaXo9Hezj6bS+ra2tyMvLwzfffIP6+vrBhumxtrY2\nHD16FDk5OQF9/rqjowOFhYXIycnx++OD3v7MvD0u9W/ChAl44IEHEBkZ6beep87OTpSVlWHPnj1e\nrW9w/fXXY8mSJSgrK/NoDImvzJs3D3fffTdqampQVlbm0XvS0tJw1113ISUlBXl5eX6OMCj9btWq\nVS5/UfHWQICYTCakpKSgubk56J7/DmU5OTnIyckJ+Hlzc3ORm5sb8POeOHECJ06cCPh5SRlxcXGY\nPn26xx94A1FRUYG3337b6/eNGjUKWVlZ+Pe//+2HqFyz2WyYNWsWNm3a5PF7oqKiMHXqVFUv6ewO\n5xEIkIsuugiPPvooDhw4gN///vdKh0NEISAmJgajRo1CRUVF0N17T09Px5AhQ1BQUBDQ2S1TUlKQ\nkpKCgoICj6f/tlqtGDlyJBobG1FQUODnCIMS5xEIBna7HcXFxaisdPmzCEkZGRmYNGkSdu3apdb/\nYER+U1tbi7179yp2/tGjRyM7Oxt79uw5b14ApQbcnTp16ryBwePGjcPEiRNd/h5qbGz0++DavpjN\nZkyfPh0OhwO7du1Ca2trwGPwiNKDA9UyWFCj0UiDwSD1er3isfiy3H///bKqqkreeOONisfCwsLi\n23LrrbfKhoYGeeeddyoei7vywAMPyMrKSrl8+XLFY+lZhg8fLj/44AP5zjvvyKSkJCVjcTtYkLcG\naFAmTZqEmTNnYvPmzThy5IjS4RCFlfHjx+NHP/oRvv76a3z66acBP/+ECRMwd+5cfPnll4r8Re2p\nyZMnY8aMGUH3e8hqtWLu3LlwOBzYtGmTkgs7ub01wESAwpJOp8OwYcMgpUR5eTk6OzuVDmnAEhMT\nYbVaUV5eHpBn700mE4YOHYrm5mZUVFT4/Xzk2mWXXYZVq1bhvffew5o1a5QOh0KX20RAE8hIiAIl\nJiYG9957L+655x5ER0crHc6gLFmyBM888wwyMjICcr4LLrgAv/rVr7Bs2bKAnI9cO3jwIFauXBnU\nU/ZS6ONgQQpLXV1dqK6uhsPhQFdXV0DOOXr0aEydOhW7d+/26SN+9fX1KC8vR1tbm8+O6U5bWxsq\nKioCOi8D9a26utqjtTmIBkXpwYFqGSzIEtgihJCRkZEyMjIyYDPY/cd//Iesq6uTt956q0+PazAY\npNlsllqtNiDfh1arlWazWRoMBsV/jiwsLD4pHCxIFAgXX3wx5s2bhy+++EKRyYaIiFzgYEEiIiKl\nxcbGIioqCpWVlYFedIuDBYmIiJT2wx/+EL/97W+RmZmpdCjnYCJA/UpKSsLVV1+Niy++WOlQiFQl\nOTkZCxcuxIQJE5QOhXygtbUVtbW16OjoUDqUcyk9OJCDBYO/XHrppXLfvn1y9erVisfCwqKmctVV\nV8m8vDz52GOPKR4Ly+CL0WiUQ4YMUWKGWbeDBfn4oAI0Gg0WL16MoUOH4sMPPwz61QiLiorwyiuv\nnDfXeCi54oorkJ2djQ8++CCoZh7rS3x8PBYvXoympiZ8+OGHsNvtAT1/WloaFi9ejOPHj2PDhg0B\nPTed6/jx43jppZeCelY/8pzdbg/4/2dP8NaAAoQQyM7Oxvz58xEbG6t0OP06efIkXn75ZWzZskXp\nUAZswoQJWLx4MYYPH650KP06Oy3pzJkzEREREfDzJyQk4JprrsHkyZMDfm461/Hjx7F27Vps375d\n6VBoEEwmE2w2W9BObsanBhQghMC4ceNgtVpx+PDhgC7hqVajRo2CzWbD4cOHceZMcF9eZrMZ48aN\nQ1tbGw4fPhzw6ZFjYmIwbtw4VFdX4+jRowE9N1E4mjFjBh544AF8/PHHWL9+vRIhcBniYCOlRF5e\nntJhqMqJEyd8Ottff6KjozFp0iTU1dXhwIED8Cbhbm5uxp49e/wYnXu1tbXYsWOHYucnCjednZ1o\naWkJvkGCTuwRIPKDzMxMPP/88zh48CAeeeQROBwOpUMiIoUYDAZYrVa0trYGev6As9gjQBRolZWV\nePvtt1FeXu62NyAxMRELFy5EfX09/vd//xft7e0BjJKIAqGtrS1ga4UMBBMB6pfBYEBMTAxaW1u5\nEI2HysrK8Le//a3fetHR0bj66qtRWlqKTz75JKwTAZPJhJiYGNTX16Opqem8/RqNBrGxsZBSora2\nlr0oRAHCpwaoXxkZGXjmmWdw3XXXKR1K2CkrK8Pzzz+PdevWBfVfDL6QnZ2N//mf/8G8efP63B8R\nEYEHH3wQd9xxB8xmc4CjI1Iv9ghQv6SUcDgcXg14CxZJSUnIzMxEcXFxUM6D0NjYiJycHKXD6FdM\nTAwyMzNRU1ODgwcPuqwXFRWFzMxMNDY24uDBg+f9Ve/uOoqOju45sRgRBYrSMwlyZsHgLwaDQSYm\nJsqoqCjFY/G2zJ8/X+7bt0/+/Oc/VzyWUC5TpkyRGzdulM8++6zbeuPHj5cbNmyQv//978+bPS0y\nMlIOGzZMms3mPt+r0+lkXFycjImJkRqNRvHvmYUljApnFqTBaWtrC/pn710pKirCP/7xD+zbt0/p\nUELamTNn8N5776G4uNhtverqanz00Uc4c+bMeb0Bra2taG1tdfnezs5OVFdX+yReIvIcHx8MEL1e\njyFDhqC9vZ0TCKmUTqfDkCFD0NnZiYaGBnaBE1GgcBniYDBq1Cg8+eSTWL58udKhkEKGDx+OX/zi\nF7jlllug07EzjoiCAxOBABFCQK/XQ6vVKh0KKUQIAZ1OxySAiIIKbw0EiMFgQGxsLOx2O2pra5UO\nhxSg1+sRGxuLzs5O1NTU8NYAEQWK21sDTARo0PR6PX7wgx8gOjoan376KSorXV5vREQUeJximPxL\np9NhwYIFSEtLw+7duwOeCGi12u8moGlqauKMdEREXmAiQIPW3t6OV199FWazGadPnw74+ePi4nDH\nHXfAbrfjlVde4VMZREReYCJAg9bV1aXosrkajQZmsxlarRZCCMXiICIKRRwjQCFPp9MhPj4eUkpU\nVVWhq6tL6ZCIiIKJ7+YREEJECCF+I4ToEEKM6GP/nUKIvUKIr4QQG4QQyX3UeUoI8Y0Q4mshxOtC\niChvYiDqrbOzE+Xl5aioqGASQETkJY8TASFEGoCtAIYCOO9heCHEUgBPAZgvpbwEwNcANogefbVC\niIcALAUwU0qZDaAdwOsDD59CRUREBIYMGYKIiAiP6huNRlgsFj5zT0TkZ970CJgB/BjAOhf7nwDw\nqpSyxvn6BQATAFwDAEIILYCVAF6UUp5db/V5AAuFEBO8DZxCy5w5c/Dyyy9j9uzZHtW/4YYbsGrV\nKowePdrPkRERqZvHiYCUMk9KWQDgvNFYQohYABMB7O1RvwFAPoCzi49nAojvWQfAUQDNPepQmIqM\njERSUhIiIyM9qm+1WhEXFwe9Xu/nyIiI1M1X/a7pzq8VvbaXA0hz/ntk7zpSSimEqOhRh8LUtm3b\n8O2336Kiovcl0re3334bGzZs8Lg+ERENjK8SAZPza1uv7e099rmq09ZjH4Wp+vp61NfXe1y/srKS\nMxQSEQWArxYdanF+NfTabuixz5M6YUGr1SIyMtLjgXFERERK8VUiUOD8mtRr+9Ae+1zVSeqxLyyk\npKRg1apVWLZsmdKhEBERueWTREBKWQtgP4CpZ7c55wcYDWCjc1MugMpedTLQ/TTCRoSRmJgYJCcn\nIyYmRulQiIiI3PJ6ZkEhxGUANgNIl1IW9di+BMAaAJlSyhohxBMAlkgpp/So8zMANwOYIaW0CyH+\nCiBeSrm4R52Qn1nQZDIhOTkZTU1NKC8vVzocIiJSN9+sPiiE0AP4HEA0AAngLSFEiZRyGQBIKT8Q\nQiQC+EIIYQdQA2Bhz2NIKf8ghLAA2CGE6ET344U3e/sdBbuWlhYcP35c6TCIiIj6xbUGiIiIwpvv\n1hogIiKi8MJEgIiISMWYCBAREakYEwEiIiIVYyKgAhqNBhoNf9RERHQ+LvYe5uLi4rBixQq0t7fj\n9ddfR2Njo9IhERFREGEi4EdJSUlISUlBUVERqqurFYkhIiICI0eORFtbG3Q6/riJiOhc7C/2o7lz\n52LNmjWYMWOGYjFUV1djzZo1ePnll9kbQERE5+GfiH505swZ7N+/f1C9AbNmzcKwYcOwfft2VFRU\neP3+9vZ2HDt2bMDnJyKi8MaZBf1MCIHBtPGaNWswZ84c3H333di5c6cPIyMiIpXwzVoDNDCDTbTe\nfvttbN++HQUFYbVSMxERBQn2CBAREYU3rjVAREREfWMiQEREpGJMBIiIiFSMiQAREZGKMREgIiJS\nMSYCREREKsZEgIiISMWYCBAREakYEwEiIiIVYyJARESkYkwEiIiIVIyJABFRCBo3bhyWLl0Km82m\ndCgU4pgIEBGFoB/84Af47//+b0yZMkXpUCjEcRliIqIQtGPHDrS1teHIkSNKh0IhjssQExERhTcu\nQ0xERER9YyJARESkYkwEiIiIVIyDBYmIQoDFYsHixYuh1+vxwQcfoK6uTumQKEywR4CIKAQYjUbM\nmTMHl112Gcxmc591oqKiYLPZXO4n6gt7BIiIQkBDQwNeeuklaLVa1NTU9Fln/vz5WL58OdasWYOt\nW7cGNkAKWUwEiEgVRo8ejbFjx2L//v0oLS1VOhyvtbe3Y//+/W7rdHZ2wm63o6urK0BRUTjgrQEi\nUoUrr7wSf/zjHzFt2jSlQ/Gbzz//HA8++CD27NmjdCgUQtgjQESqsH//frzyyivIz89XOhS/aWlp\nQUtLi9JhUIjhzIJE5LXo6Gjo9XrU1dWho6PjvP0mkwnR0dGor69Hc3OzAhGGF6vVCrPZjLq6Otjt\ndqXDodDDmQWJyHc0Gg1WrFiBxx9/HKmpqX3WmTVrFl588UVcfvnlAY4uPF1zzTX43e9+h4kTJyod\nCoUh3hogogFx15sopXS732KxIDMzE3a7Hbm5uejs7PRHiGEl2HpvKXzw1gAReW2wtwZGjRqF3/zm\nNzhz5gxWrlyJxsbGQIQdsnhrgAbJ7a0B9ggQkdf6m9Wuv0Fr9fX1+OSTT9DY2NhnIkHnamxsHFSy\npNfrcfXVVyM2NhYff/wxzpzh31r0PSYCRD4SGRmJqKgoNDQ0oLW1VelwvCaEQFRUFDQaDRobG/3a\nXV9VVYV169b57fh9MZlMsFgsaGhoCLm/qs1mM0wmExoaGtDW1ub1+3U6HebNm4e0tDTk5OQwEaBz\ncLAgkY9ceumleOmll3DJJZcoHcqAGAwG3HbbbXjkkUcwbNgwpcPxuVCeR2CwgwXb29uxbt06/Pa3\nv8Xp06eJSHiIAAAT0UlEQVR9HB2FOvYIEPmIVquF0WiEVqv1qH5GRgbi4uKQl5eH2tra8/YPGzYM\n48aNQ0FBAQoLC30d7nmEENDr9YiIiIAQwu/nCzSdTgej0QiNJvT+/tHpdDAYDAOOvaurC/v27fNx\nVBQ2zo7uDZYCIAGAZGEJtWK1WmVqaqq0WCwe1V+1apX8/PPPZVZWVp/7Fy5cKHNzc+Vdd90VkPiF\nEDI+Pl4OHTpU6vV6xdvT1yUqKkrabDZpNpsVj8XbMmTIEDl8+HBpMpkUj4UlJEuCu89d9ggQ9WPa\ntGmYM2cO/v3vf+PQoUMu63k7oGvv3r2ora1FZWXfg3mLiorw7rvv4vDhw17HPBBSSlRVVQXkXEpo\naGhAQ0OD0mEMSH19Perr65UOg8IUEwGifkyZMgUPPPAAiouL3SYC3tqwYYPb/bm5ucjNzfXZ+YjO\nMhqNMJvNaGpqGtDgQwovoXezjCjANm3ahJ/+9Kf4+uuvlQ6FyCcuv/xyrF27FjNnzlQ6FAoC7BEg\n6kd+fn5YL1RD4W3EiBEYNWoUjh49irKyMgDdT4jExMTAYDAoHB0FA/YIEBGFsSuuuAJr1qzB9OnT\nv9u2ZcsW3HPPPdixY4eCkVGwYI8AEZ1j+PDhmD9/PoqKirBp0yalw/GZrKwszJo1C5999hmOHDmi\ndDhuWSwWLFiwAFJKfP7552hqanJbPzU1FfPnz8exY8ewdevWc/YdP34cH374IYqKir7bxsGH1BN7\nBIjoHDabDXfeeSeuvPJKpUPxqenTp+Phhx/G+PHjlQ6lX1arFcuWLcOyZctgtVr7rZ+eno577rkH\nc+fOPW/f9u3b8cQTT3AeAXKJiw4R0TkSEhKQlZWF8vLysPrwGDt2LCZMmIA9e/ac89dxMIqMjERW\nVhYAYM+ePf1OWT106FBkZWWhpKQEBw4cCESIFFrcLjrERICIiCi8uU0EeGuAiIhIxZgIkFciIiJw\n9dVXY9myZYiNjVU6HCLykczMTPzsZz/DRRddpHQoFGB8aoC8YjAYsGzZMgwbNgwHDhxATU2N0iH5\nnU6nQ0REBNrb2/26NC+RkqZNm4ann34a9fX1OHjwoNLhUABxjAB5Ra/XIzs7GyaTCbt27QrZudu9\ncdVVV+G6667Dq6++yueuKWyNGTMG06ZNw9dff41jx44pHQ75ltsxAuwRIK90dHTgq6++UjqMgBoy\nZAjS0tJgsViUDoXIbziDpnqxR4CoH3FxcUhMTERZWRnq6uqUDoeIyFt8fJCIgsfs2bORmpqKLVu2\n4NSpUwE9d0REBC677DJYLBZs2bIFtbW1AT3/QGVlZSEzMxNbtmxBQUGB0uFQ6OGtASIaOJ1OByEE\nOjs74e4PB0/rXXvttbjqqqtQXFwc8ETAYDDgpptuwvDhw3H48OGQSQTmz5+Pu+66C5WVlUwEyOeY\nCBCRSzqdDrfccguGDh2K1157DaWlpS7rLl++HBdccAFee+01tx9W7777LnJychS5H2232/H3v/8d\nJpPpu5X4QsHHH3+MEydOcNZA8gsmAkTkkhACycnJGDFiBIxGo9u6w4YNQ3p6OiIjI93W2717N3bv\n3u3LMD3W0dGBL7/8UpFzD8aBAwcUTQLS0tIQFRWFgoKCfhdAotDDMQJE5JIQAjabDQaDAaWlpW7n\nvB8+fDjMZjNKSkrQ0tISwCjJ3x577DFMmjQJzzzzDHJzc5UOh7zHMQJEgSKEwJw5c5CQkIAvv/wS\nlZUu/++FBCkliouLParr7rYBhbbCwkJoNBo0NjYqHQr5ARMB8poQAkIIOBwOpUMJOkII3HjjjZg8\neTIKCwtDPhEIR6Fy/QohAMDtwMtAeeutt5QOgfyItwbIK0ajEStWrMCQIUOwfv16lJeXKx1SUBFC\n4JJLLkF8fDy2b9+OqqoqpUOiHoxGI3784x8jJiYG69evD9oBg9HR0VixYgWklFi/fj3nr6DB4q0B\n8h2tVou0tDQkJCT0O3gsHFksFowYMQKNjY0oKSk5b7+UEtu3b/foWDabDVarFcXFxRyAFSAajQZp\naWlISkoK6us3IiICo0aNgpQSERERSodDYY49AuQVrVaL1NRU6PV6FBUVwW63Kx1SQE2cOBGPPfYY\n9uzZg+eff35Qx3r44YeRlZWF5557jo+FBUioXL96vR6pqakAgKKiInR0dCgcEYU49giQ73R1dal6\nQpOmpiYcOXLEJwPjSktLYbFY0Nzc7IPI/CsxMRGzZ89GeXl5SC+8FCrXb0dHB44fP650GKQS7BEg\n8pJGo4GUctCDuEJl0BoATJ8+HS+//DK++uor3HfffUqHQ0Te8V2PgBAiAsB/Avg5gFFSyuIe+14F\ncCGAnn1th6WU9/Y6xlMAFgHoBJAP4F4pZfivZUthw1cf3L5IJgKlsLAQzz77bNAOrgskm82GG2+8\nEadOncKbb76Jrq4uJCcn48Ybb0RlZSXeeOMNduVTSPE4ERBCpAF4A90f3to+qkgAN/RMDvo4xkMA\nlgLIllK2CSH+BuB1dCcGRBSkKioq8M477ygdRlAwm8246KKLEBkZ+d0jfme3lZSUQKPRKBwhkXc8\nvjUghBgPoBWADcAWAGm9egTWAVglpSxy8X4tgDIAv5RS/sW5LQNAHoBMKeUh5zbeGiCioGU2m5Ge\nno6WlhYUFhZCSgmTyYT09HTY7XYUFhaGxO0eUhXf3BqQUuYBgBBixAADyQQQD2Bvj21HATQDmAfg\n0ACPS0TkUxqNBjNmzEBMTAxycnJQXV393b7m5mYcOnTur6uWlhbk5eUFOkwin/B1H9ZjQogtQojt\nQoi1QojEHvtGOr9WnN0gu7sjKgCk+TgOIqIB02g0uPnmm7Fy5UqMGDHQv32IQoMvHx/MB3BSSnmX\nEEIP4E8AdgkhLpJSNgMwOeu19XpfW499RER+N2rUKNxwww04evQo3n///fP2OxwOvPPOO4iLi+Ma\nChT2fJYISCmf6/HvDufAwFoAywH8FcDZ5cgMvd5q6LGPglh8fDzi4+NRVlaG+vp6pcMhGrDo6GhM\nnz4dXV1dfe53OBzYtGlTgKMiUobfJhSSUjYKISrx/S2Bs7N4JAE43aNqUo99FMRGjx4Ni8WC4cOH\nY+PGjUqHQzRghw4dwqpVqziHPxF8OEZACPFCr9cGAHEAzj5ZkAugEsDUHnUyAJgB8FMlBMTExOD6\n66/nvPgAxo4dixtvvBFWq1XpUGgA2trasG/fvpCYZZDI785OauJpAXAZAAeA1F7b7QCm9Hj9NIBy\nAHE9tv0MwH4ARufrvwL4sNdxEtA9JwFLkBW9Xi9NJpPU6XSKx6J0ufvuu2VNTY2cNWuW4rGwsLCw\n9FMS3H2uezOhkB7A5wCinQd+SwhRIqVc5qzyMIA/CCE60T347wyAy6WU3z13I6X8gxDCAmCHs14+\ngJs9jYGU1dHRwRnTnHbt2oXVq1eza7mXlJQUXH/99SgvL8c///nPPu/BL1q0CBkZGXj33Xdx4sQJ\nl8e6+uqrMXnyZPzzn//E0aNH/Rk2kap5M49AB4DL3exfC2CtB8d5Gt29BUQha//+/di/f7/SYQRM\nZGQkkpKS0NraioqKCpf1oqKiMGvWLBw7dgzvv/9+n4nAuHHjMGfOHGzevNltInDhhRdi7ty52LZt\nGxMBIj/iokNE1K/x48fjoYcewsGDB/HHP/7RZT2LxYIxY8agubkZ+fn5fa6lcMEFFyA2Nhb5+flu\ne1TS09ORlJSE/Px81NTU+OT7IFIpLkNMRIPT3t6OqqoqNDS4Xx+sqakJ+/btc1vH0+V1CwsLUVhY\n6HGMRDQw7BEgon5ptVpERkaiq6sLra2tSodDRN5hjwARDU5XVxcfGyUKU1wvkygExcbGIjU1FSYT\nZ+cmosFhIkAUghYtWoRnn30WEyZMUDoUIgpxTATIY2PGjMGiRYtgs9mUDiWsWSwWzJkzB9nZ2dDr\n9X3WsdvtaGxsRGdnZ4CjI6Kw4+3Mgv4u4MyCQVvuu+8+WVhYKJcuXap4LOFcRo0aJT/88EP5l7/8\nRQ4ZMqTPOiaTScbExEi9Xq94vCwsLEFffDOzINHevXuxdu1aTu4yQMOGDcOiRYtQWVmJjz76yOVf\n87W1tXj33XfR3NyMtrbeq3Z3a2lpQUsLF+0kosHj44NEAySEQGxsLHQ6HWpqavqdfjkjIwPPPfcc\nTpw4gcceewzt7e0BipSIVM7t44NMBIgGSKfT4f7774fNZsOaNWv6nfzGarViwoQJaGpqQl5eHhwO\nR4AiJSKV4zwCRP7S2dmJjo6OPqfS7a2xsRE5OTkBiIqIyHPsESAaICEETCYTDAYDGhoaOIKfiIIV\newSI/EFKiebmZjQ3NysdChHRgHEeASIiIhVjIkBERKRiTASIiIhUjIkAERGRijERICIiUjEmAjQg\ner0esbGxsFgsSodCRESDwESABmTMmDFYvXo1fvSjHykdChERDQITARoQjUYDo9HocpncQBkzZgzm\nzJmD+Ph4ReMgIgpVnFmQBsRoNCIuLg6tra2oqalRLI5f/OIXuOaaa/DEE09g27ZtisVBRBTEOLMg\n+Z7dbsepU6eUDgMHDhyAw+FAWVmZ0qEQEYUk9giQ4kwmE/R6PZqbmzlfPxGR77ntEeAYAVLc8uXL\n8eSTT2LkyJFKh0JEpDpMBEhxZrMZUVFR0Ol4p4qIKNB4a4AUl5CQAKPRiMrKStjtdqXDISIKNxws\nSMGtstLl9UlERH7GWwNEREQqxkSAiIhIxXhrgAZNo9FgzJgxMBqNOHbsGJqbm5UOiYiIPMQeARq0\niIgI3HrrrXj00UeRnJysdDhEROQF9gjQoDkcDuTm5uL06dNobGxUOhwiIvICHx8kIiIKb5xZkIiI\niPrGRICIiEjFmAgQERGpGBMBIiIiFWMiQEREpGJMBIiIiFSMiQAREZGKMREgIiJSMSYCREREKhaM\niYBQOgAiIqIw4vZzNRgTgVilAyAiIgojbj9XgzERICIiogBhIkBERKRiTASIiIhULBiXIdYCGN1r\ncw2A4AqUiIgo+AicPybgmJSyy+Ubgi0RICIiosDhrQEiIiIVYyJARESkYjqlA3BHCLEEwGMA7AAc\nAO6RUh5WNqrgJoRYBWARgLoem2uklNf1qHMngJ+gu13rANwhpTwdyDiDmRAiAsB/Avg5gFFSyuJe\n+/ttPyHEU+j+OXQCyAdwr5SyIQDhBxV3bSmEeBXAhehux7MOSynv7XUMtqWTEGIZgP8AoAUQBeAk\ngEeklEU96vD69EB/bamq61NKGZQFwDQA9ej+5QEAKwCUALAoHVswFwC/AjDHzf6lAE4BiHW+fhLA\nPjjHi6i9AEgDsBPAq+hOPkd4234AHgJwAIDB+fpvAD5S+nsLwrZc13tbH8dgW57bHm0ArnT+WwB4\nDcARABHObbw+B9+Weuc21VyfigfgpoHfB/CPHq8FgDIA9ykdWzAXZyJwqZv9+wD8usfrKADtAH6o\ndOzBUACMBzASwKUuPrzcth+6/7o4A+AnPepkOI81QenvL8jach2AVDfvZ1ue3yZv93o9xdke2c7X\nvD5915aquT6DeYzAXAB7z76Q3a38DYB5ikUU4oQQsQAm4tx2bUB3dxbbFYCUMk9KWYA+5ub2sP0y\nAcT3rAPgKIBmqKyN3bWlh9iWvUgpb+i1qc351cDr0zvu2tLDQ4RNWwZlIiCEiEN3JlvRa1cFgPTA\nRxRybhNCbBFCfCWEeFUIMdK5/Wzb9W7XcrBdPeGu/dKc/x7Zu44zia3oUYe+95jzWt0uhFgrhEjs\nsY9t2b8Z6L4VsAO8PgerZ1uepYrrMygTAQAm59e2XtvbeuyjvhWju3vwCgCzARQC+EYIkQzX7doO\nIDJgEYYud+1n6qcOr93z5QP4Ukp5Obp7AA0AdgkhzM79bEs3hBAGAA+j+3ZpF3h9DlgfbQmo6PoM\n1kSgxfm1dxeNAd3dLuSClHKdlPIFKaXDmZ0+je5Rr/fg+7brq11bQP1xd122eFGHAEgpn5NSvun8\ndwe6B16NALDcWYVt6d6fAbwlpfzI+ZrX58D1bktVXZ9BmQhIKavR/cRAUq9dQwEUBD6i0CWldKD7\nsZiR+L7t2K4D40n7uaqTBLaxW1LKRgCV+L7LlW3pghDiNwCapJS/6rGZ1+cAuGjL84Tz9RmUiYDT\nZgBTz74QQggAkwFsVCyiECCEeKGPzSkAiqWUdQD249x2jUL32g5s135IKWvRf/vlovuXRc86GQDM\nYBufo/e16uyejUP37S2AbdknIcRKAClSyvucrycLISbz+vSeq7Z0/ls916fSjy24KgCy0D0Zxtl5\nBH6M7h+AWenYgrmgOxNd2OP17ei+JTDG+XoJgFJ8/5zxEwC+UTruYCsALkP3Y0Cpvbb3234Afobu\nX8hG5+u/AvhQ6e8pCNvSDmBKj9dPo3tgWxzb0mVb3gXgIIBsdH8ATQWwCsAtzv28Pn3Xlqq5PoN2\nZkEp5R4hxK0A3hJCtALoArBASskxAu79EsCDQoiHAESg+2KeJ6XMBwAp5QfOka9fCCHs6F7ZcaFi\n0QYZIYQewOcAotG94uVbQogSKeUywLP2k1L+QQhhAbBDCHF2trGbA/l9BIP+2hLdg7P+4GwjE7qf\nyb5cdt8aBMC27EkIYQWwFt2PY+b02CUB/D+A16enPGlLqOj65OqDREREKhbMYwSIiIjIz5gIEBER\nqRgTASIiIhVjIkBERKRiTASIiIhUjIkAERGRijERICIiUjEmAkRERCrGRICIiEjFmAgQERGpGBMB\nIiIiFWMiQEREpGL/H2jpzlqIoc6YAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1042eadd0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "#From http://stackoverflow.com/questions/36372068/how-to-get-x-y-coordinates-using-python\n",
    "#\n",
    "#Took picture of Leela's constellation from google\n",
    "\n",
    "import cv2\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "filename = 'Leelas_constellation.png'\n",
    "img = cv2.imread(filename)\n",
    "# plt.imshow(img)\n",
    "\n",
    "img=cv2.cvtColor(img,cv2.COLOR_BGR2GRAY) #converting to grayscale\n",
    "img=img.astype(np.uint8)\n",
    "cp_img = img.copy()\n",
    "img[cp_img == 255] = 0\n",
    "img[cp_img == 0] = 1\n",
    "plt.imshow(img, cmap='gray')\n",
    "\n",
    "#get all non zero values\n",
    "img = np.flipud(img)\n",
    "coord=cv2.findNonZero(img)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "Then create REBOUND simulation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "sim = rebound.Simulation()\n",
    "sim.integrator = \"ias15\" \n",
    "sim.collision = \"direct\"\n",
    "sim.collision_resolve = \"merge\"\n",
    "\n",
    "mass = 1.\n",
    "\n",
    "for i in range(len(coord)):\n",
    "    sim.add(m=mass, x=coord[i][0][0], y=coord[i][0][1], z=0)    \n",
    "    \n",
    "sim.move_to_com()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "\n",
       "<script id=\"orbit_shader-vs\" type=\"x-shader/x-vertex\">\n",
       "    uniform vec3 focus;\n",
       "    uniform vec3 aef;\n",
       "    uniform vec3 omegaOmegainc;\n",
       "    attribute float lintwopi;\n",
       "    varying float lin;\n",
       "    uniform mat4 mvp;\n",
       "    const float M_PI = 3.14159265359;\n",
       "    void main() {\n",
       "       float a = aef.x;\n",
       "       float e = aef.y;\n",
       "       float f = aef.z+lintwopi;\n",
       "       lin = lintwopi/(M_PI*2.);\n",
       "       if (e>1.){\n",
       "           float theta_max = acos(-1./e);\n",
       "           f = 0.0001-theta_max+1.9998*lin*theta_max;\n",
       "           lin = sqrt(min(0.5,lin));\n",
       "       }\n",
       "       float omega = omegaOmegainc.x;\n",
       "       float Omega = omegaOmegainc.y;\n",
       "       float inc = omegaOmegainc.z;\n",
       "       float r = a*(1.-e*e)/(1. + e*cos(f));\n",
       "       float cO = cos(Omega);\n",
       "       float sO = sin(Omega);\n",
       "       float co = cos(omega);\n",
       "       float so = sin(omega);\n",
       "       float cf = cos(f);\n",
       "       float sf = sin(f);\n",
       "       float ci = cos(inc);\n",
       "       float si = sin(inc);\n",
       "       vec3 pos = vec3(r*(cO*(co*cf-so*sf) - sO*(so*cf+co*sf)*ci),r*(sO*(co*cf-so*sf) + cO*(so*cf+co*sf)*ci),+ r*(so*cf+co*sf)*si);\n",
       "       gl_Position = mvp*(vec4(focus+pos, 1.0));\n",
       "    }\n",
       "</script>\n",
       "<script id=\"orbit_shader-fs\" type=\"x-shader/x-fragment\">\n",
       "    precision mediump float;\n",
       "    varying float lin;\n",
       "    void main() {\n",
       "      float fog = max(max(0.,-1.+2.*gl_FragCoord.z),max(0.,1.-2.*gl_FragCoord.z));\n",
       "      gl_FragColor = vec4(1.,1.,1.,sqrt(lin)*(1.-fog));\n",
       "    }\n",
       "</script>\n",
       "<script id=\"point_shader-vs\" type=\"x-shader/x-vertex\">\n",
       "    attribute vec3 vp;\n",
       "    uniform mat4 mvp;\n",
       "    //uniform vec4 vc;\n",
       "    //varying vec4 color;\n",
       "    void main() {\n",
       "      gl_PointSize = 15.0;\n",
       "      gl_Position = mvp*vec4(vp, 1.0);\n",
       "      //color = vc;\n",
       "    }\n",
       "</script>\n",
       "<script id=\"point_shader-fs\" type=\"x-shader/x-fragment\">\n",
       "    precision mediump float;\n",
       "    //varying vec4 color;\n",
       "    void main() {\n",
       "      vec2 rel = gl_PointCoord.st;\n",
       "      rel.s -=0.5;\n",
       "      rel.t -=0.5;\n",
       "      if (length(rel)>0.25){\n",
       "         gl_FragColor = vec4(0.,0.,0.,0.); \n",
       "      }else{\n",
       "         vec4 cmod = vec4(1.,1.,1.,1.);\n",
       "         float fog = max(max(0.,-1.+2.*gl_FragCoord.z),max(0.,1.-2.*gl_FragCoord.z));\n",
       "         cmod.a*= (1.-fog)*min(1.,1.-4.*(length(rel)/0.25-0.75));\n",
       "         gl_FragColor = cmod;\n",
       "      }\n",
       "    }\n",
       "</script>\n",
       "\n",
       "<script>\n",
       "function compileShader(glr, shaderSource, shaderType) {\n",
       "  // Create the shader object\n",
       "  var shader = glr.createShader(shaderType);\n",
       " \n",
       "  // Set the shader source code.\n",
       "  glr.shaderSource(shader, shaderSource);\n",
       " \n",
       "  // Compile the shader\n",
       "  glr.compileShader(shader);\n",
       " \n",
       "  // Check if it compiled\n",
       "  var success = glr.getShaderParameter(shader, glr.COMPILE_STATUS);\n",
       "  if (!success) {\n",
       "    // Something went wrong during compilation; get the error\n",
       "    throw \"could not compile shader:\" + glr.getShaderInfoLog(shader);\n",
       "  }\n",
       " \n",
       "  return shader;\n",
       "}\n",
       "function createShaderFromScript(glr, scriptId, opt_shaderType) {\n",
       "  // look up the script tag by id.\n",
       "  var shaderScript = document.getElementById(scriptId);\n",
       "  if (!shaderScript) {\n",
       "    throw(\"*** Error: unknown script element\" + scriptId);\n",
       "  }\n",
       " \n",
       "  // extract the contents of the script tag.\n",
       "  var shaderSource = shaderScript.text;\n",
       " \n",
       "  // If we didn't pass in a type, use the 'type' from\n",
       "  // the script tag.\n",
       "  if (!opt_shaderType) {\n",
       "    if (shaderScript.type == \"x-shader/x-vertex\") {\n",
       "      opt_shaderType = glr.VERTEX_SHADER;\n",
       "    } else if (shaderScript.type == \"x-shader/x-fragment\") {\n",
       "      opt_shaderType = glr.FRAGMENT_SHADER;\n",
       "    } else if (!opt_shaderType) {\n",
       "      throw(\"*** Error: shader type not set\");\n",
       "    }\n",
       "  }\n",
       " \n",
       "  return compileShader(glr, shaderSource, opt_shaderType);\n",
       "};\n",
       "function createProgramFromScripts( glr, vertexShaderId, fragmentShaderId) {\n",
       "  var vertexShader = createShaderFromScript(glr, vertexShaderId, glr.VERTEX_SHADER);\n",
       "  var fragmentShader = createShaderFromScript(glr, fragmentShaderId, glr.FRAGMENT_SHADER);\n",
       "  var program = glr.createProgram();\n",
       " \n",
       "  // attach the shaders.\n",
       "  glr.attachShader(program, vertexShader);\n",
       "  glr.attachShader(program, fragmentShader);\n",
       " \n",
       "  // link the program.\n",
       "  glr.linkProgram(program);\n",
       " \n",
       "  // Check if it linked.\n",
       "  var success = glr.getProgramParameter(program, glr.LINK_STATUS);\n",
       "  if (!success) {\n",
       "      // something went wrong with the link\n",
       "      throw (\"program filed to link:\" + glr.getProgramInfoLog (program));\n",
       "  }\n",
       " \n",
       "  return program;\n",
       "}\n",
       "function quat2mat(A,mat){\n",
       "    var xx = A.x*A.x; var xy = A.x*A.y; var xz = A.x*A.z;\n",
       "    var xw = A.x*A.w; var yy = A.y*A.y; var yz = A.y*A.z;\n",
       "    var yw = A.y*A.w; var zz = A.z*A.z; var zw = A.z*A.w;\n",
       "    mat[0] = 1.-2.*(yy+zz);\n",
       "    mat[1] =    2.*(xy-zw);\n",
       "    mat[2] =    2.*(xz+yw);\n",
       "    mat[4] =    2.*(xy+zw);\n",
       "    mat[5] = 1.-2.*(xx+zz);\n",
       "    mat[6] =    2.*(yz-xw);\n",
       "    mat[8] =    2.*(xz-yw);\n",
       "    mat[9] =    2.*(yz+xw);\n",
       "    mat[10]= 1.-2.*(xx+yy);\n",
       "    mat[3] = mat[7] = mat[11] = mat[12] = mat[13] = mat[14] = 0.; mat[15]= 1.;\n",
       "}\n",
       "function multvec(A, B, vecr){\n",
       "    var mat = [0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.];\n",
       "    quat2mat(A,mat);\n",
       "    vecr[0] = mat[0]*B[0] + mat[1]*B[1] + mat[2]*B[2];\n",
       "    vecr[1] = mat[4]*B[0] + mat[5]*B[1] + mat[6]*B[2];\n",
       "    vecr[2] = mat[8]*B[0] + mat[9]*B[1] + mat[10]*B[2];\n",
       "}\n",
       "function mattransp(mat){\n",
       "    var matt = [\n",
       "        mat[0], mat[4], mat[8], mat[12],\n",
       "        mat[1], mat[5], mat[9], mat[13],\n",
       "        mat[2], mat[6], mat[10], mat[14],\n",
       "        mat[3], mat[7], mat[11], mat[15]];\n",
       "    return matt;\n",
       "}\n",
       "function conjugate(quat){\n",
       "    var cquat = {x:-quat.x, y:-quat.y, z:-quat.z, w:quat.w};\n",
       "    return cquat;\n",
       "}\n",
       "function mult(A, B){\n",
       "    var mquat = {   x: A.w*B.x + A.x*B.w + A.y*B.z - A.z*B.y,\n",
       "                    y: A.w*B.y - A.x*B.z + A.y*B.w + A.z*B.x,\n",
       "                    z: A.w*B.z + A.x*B.y - A.y*B.x + A.z*B.w,\n",
       "                    w: A.w*B.w - A.x*B.x - A.y*B.y - A.z*B.z};\n",
       "    return mquat;\n",
       "}\n",
       "\n",
       "function normalize(quat){\n",
       "    var L = Math.sqrt(quat.x*quat.x + quat.y*quat.y + quat.z*quat.z + quat.w*quat.w);\n",
       "    var nquat = {x:quat.x/L, y:quat.y/L, z:quat.z/L, w:quat.w/L};\n",
       "    return nquat;\n",
       "}\n",
       "function matortho(mat, l, r, b, t, n, f){\n",
       "    mat[0] = 2./(r-l); mat[1] = 0.; mat[2] = 0.; mat[3] = -(r+l)/(r-l);\n",
       "    mat[4] = 0.; mat[5] = 2./(t-b); mat[6] = 0.; mat[7] = -(t+b)/(t-b);\n",
       "    mat[8] = 0.; mat[9] = 0.; mat[10] = -2./(f-n); mat[11] = -(f+n)/(f-n);\n",
       "    mat[12] = 0.; mat[13] = 0.; mat[14] = 0.; mat[15] = 1.;\n",
       "}\n",
       "function matmult(A,B,C){\n",
       "    for(i=0;i<4;i++){\n",
       "    for(j=0;j<4;j++){\n",
       "        C[i+4*j] = 0.;\n",
       "    for(k=0;k<4;k++){\n",
       "        C[i+4*j] += A[k+4*j]*B[i+4*k];\n",
       "    }}}\n",
       "}\n",
       "function startGL(reboundView) {\n",
       "    var canvas = document.getElementById(\"reboundcanvas-\"+reboundView.id);\n",
       "    if (!canvas){\n",
       "        reboundView.startCount = reboundView.startCount+1;\n",
       "        if (reboundView.startCount>1000){\n",
       "            console.log(\"Cannot find element.\");\n",
       "        }else{\n",
       "            setTimeout(function(){ startGL(reboundView); }, 10);\n",
       "        }\n",
       "        return;\n",
       "    }\n",
       "    var rect = canvas.getBoundingClientRect()\n",
       "    reboundView.ratio = rect.width/rect.height;\n",
       "    reboundView.view = normalize({x:reboundView.orientation[0], y:reboundView.orientation[1], z:reboundView.orientation[2], w:reboundView.orientation[3]});\n",
       "\n",
       "    canvas.addEventListener('mousedown', function() {\n",
       "        reboundView.mouseDown=1;\n",
       "        }, false);\n",
       "    canvas.addEventListener('mouseup', function() {\n",
       "        reboundView.mouseDown=0;\n",
       "        }, false);\n",
       "    canvas.addEventListener('mouseleave', function() {\n",
       "        reboundView.mouseDown=0;\n",
       "        }, false);\n",
       "\n",
       "    canvas.addEventListener('mousemove', function(evt) {\n",
       "        var rect = canvas.getBoundingClientRect()\n",
       "        if (reboundView.mouseDown==1){\n",
       "            reboundView.mouseDown = 2;\n",
       "            reboundView.mouse_x = evt.clientX-rect.left;\n",
       "            reboundView.mouse_y = evt.clientY-rect.top;\n",
       "            return;\n",
       "        }else if (reboundView.mouseDown==2){\n",
       "            var width = rect.width;\n",
       "            var height = rect.height;\n",
       "            var dx = 3.*(evt.clientX-rect.left-reboundView.mouse_x)/width;\n",
       "            var dy = 3.*(evt.clientY-rect.top-reboundView.mouse_y)/height;\n",
       "            reboundView.mouse_x = evt.clientX-rect.left;\n",
       "            reboundView.mouse_y = evt.clientY-rect.top;\n",
       "            if (evt.shiftKey){\n",
       "                reboundView.scale *= (1.+dx+dy);\n",
       "            }else{\n",
       "                var inv = conjugate(reboundView.view);\n",
       "                var up = [0.,1.,0.];\n",
       "                var right = [1.,0.,0.];\n",
       "                var inv_up = [0.,0.,0.];\n",
       "                var inv_right = [0.,0.,0.];\n",
       "                multvec(inv, right, inv_right);\n",
       "                multvec(inv, up, inv_up);\n",
       "                \n",
       "                var sin_dy = Math.sin(dy);\n",
       "                var rot_dy = {x:inv_right[0]*sin_dy, y:inv_right[1]*sin_dy, z:inv_right[2]*sin_dy, w:Math.cos(dy)};\n",
       "                reboundView.view = mult(reboundView.view, normalize(rot_dy));\n",
       "                \n",
       "                var sin_dx = Math.sin(dx);\n",
       "                var rot_dx = {x:inv_up[0]*sin_dx, y:inv_up[1]*sin_dx, z:inv_up[2]*sin_dx, w:Math.cos(dx)};\n",
       "                reboundView.view = normalize(mult(reboundView.view, normalize(rot_dx)));\n",
       "            }\n",
       "\n",
       "            drawGL(reboundView);\n",
       "        }\n",
       "\n",
       "\n",
       "        }, false);\n",
       "\n",
       "    reboundView.gl = canvas.getContext(\"webgl\")||canvas.getContext(\"experimental-webgl\");\n",
       "    if (!reboundView.gl) {\n",
       "      alert(\"Unable to initialize WebGL. Your browser may not support it.\");\n",
       "      return;\n",
       "    }\n",
       "    var gl = reboundView.gl\n",
       "    gl.enable(gl.BLEND);\n",
       "    gl.blendFunc(gl.SRC_ALPHA, gl.ONE_MINUS_SRC_ALPHA);\n",
       "    \n",
       "    reboundView.orbit_shader_program = createProgramFromScripts(gl,\"orbit_shader-vs\",\"orbit_shader-fs\");\n",
       "    reboundView.point_shader_program = createProgramFromScripts(gl,\"point_shader-vs\",\"point_shader-fs\");\n",
       "   \n",
       "    var lintwopi = new Float32Array(500);\n",
       "    for(i=0;i<500;i++){\n",
       "        lintwopi[i] = 2.*Math.PI/500.*i;\n",
       "    }\n",
       "    reboundView.orbit_lintwopi_buffer = gl.createBuffer();\n",
       "    gl.bindBuffer(gl.ARRAY_BUFFER, reboundView.orbit_lintwopi_buffer);\n",
       "    gl.bufferData(gl.ARRAY_BUFFER, 4*500, gl.STATIC_DRAW);\n",
       "    gl.bufferSubData(gl.ARRAY_BUFFER, 0, lintwopi)\n",
       "    reboundView.orbit_shader_mvp_location = gl.getUniformLocation(reboundView.orbit_shader_program,\"mvp\");\n",
       "    reboundView.orbit_shader_focus_location = gl.getUniformLocation(reboundView.orbit_shader_program,\"focus\");\n",
       "    reboundView.orbit_shader_aef_location = gl.getUniformLocation(reboundView.orbit_shader_program,\"aef\");\n",
       "    reboundView.orbit_shader_omegaOmegainc_location = gl.getUniformLocation(reboundView.orbit_shader_program,\"omegaOmegainc\");\n",
       "    \n",
       "    reboundView.particle_data_buffer = gl.createBuffer();\n",
       "    gl.useProgram(reboundView.point_shader_program);\n",
       "    reboundView.point_shader_mvp_location = gl.getUniformLocation(reboundView.point_shader_program,\"mvp\");\n",
       "    \n",
       "    updateRenderData(reboundView);\n",
       "    gl.clearColor(0.0, 0.0, 0.0, 1.0);\n",
       "    gl.clear(gl.COLOR_BUFFER_BIT);\n",
       "    drawGL(reboundView);\n",
       "}\n",
       "function updateRenderData(reboundView){\n",
       "    var previousN = reboundView.N;\n",
       "    reboundView.N = reboundView.model.get(\"N\");\n",
       "    reboundView.t = reboundView.model.get(\"t\");\n",
       "    reboundView.particle_data = reboundView.model.get('particle_data');\n",
       "    if (reboundView.orbits){\n",
       "        reboundView.orbit_data = reboundView.model.get('orbit_data');\n",
       "    }\n",
       "    var gl = reboundView.gl\n",
       "    if (reboundView.N>0){\n",
       "        gl.bindBuffer(gl.ARRAY_BUFFER, reboundView.particle_data_buffer);\n",
       "        gl.bufferData(gl.ARRAY_BUFFER, reboundView.N*7*4, gl.DYNAMIC_DRAW);\n",
       "        gl.bufferSubData(gl.ARRAY_BUFFER, 0, reboundView.particle_data)\n",
       "    }\n",
       "}\n",
       "function drawGL(reboundView) {\n",
       "    if (!reboundView.gl){\n",
       "        return;\n",
       "    }\n",
       "    // Cleanup\n",
       "    var gl = reboundView.gl\n",
       "    gl.clearColor(0.0, 0.0, 0.0, 1.0);\n",
       "    gl.clear(gl.COLOR_BUFFER_BIT);\n",
       "    \n",
       "    // Draw\n",
       "    gl.useProgram(reboundView.point_shader_program);\n",
       "    gl.bindBuffer(gl.ARRAY_BUFFER, reboundView.particle_data_buffer);\n",
       "    var pvp = gl.getAttribLocation(reboundView.point_shader_program,\"vp\");\n",
       "    gl.enableVertexAttribArray(pvp);\n",
       "    gl.vertexAttribPointer(pvp, 3, gl.FLOAT, 0, 4*7,0); // 4 = size of float\n",
       "    var projection = [0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.];\n",
       "    if (reboundView.ratio>=1.){\n",
       "        matortho(projection, \n",
       "                -1.6*reboundView.scale, 1.6*reboundView.scale,\n",
       "                -1.6/reboundView.ratio*reboundView.scale, 1.6/reboundView.ratio*reboundView.scale,\n",
       "                -2.5*reboundView.scale, 2.5*reboundView.scale);\n",
       "    }else{\n",
       "        matortho(projection, \n",
       "                -1.6*reboundView.ratio*reboundView.scale, 1.6*reboundView.ratio*reboundView.scale,\n",
       "                -1.6*reboundView.scale, 1.6*reboundView.scale,\n",
       "                -2.5*reboundView.scale, 2.5*reboundView.scale);\n",
       "    }\n",
       "    var view = [0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.];\n",
       "    quat2mat(reboundView.view,view);\n",
       "    var mvp = [0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.];\n",
       "    matmult(projection,view,mvp);\n",
       "    gl.uniformMatrix4fv(reboundView.point_shader_mvp_location,false,mattransp(mvp));\n",
       "    gl.drawArrays(gl.POINTS,0,reboundView.N);\n",
       "   \n",
       "    if (reboundView.orbits){\n",
       "        gl.useProgram(reboundView.orbit_shader_program);\n",
       "        gl.bindBuffer(gl.ARRAY_BUFFER, reboundView.orbit_lintwopi_buffer);\n",
       "        var ltp = gl.getAttribLocation(reboundView.orbit_shader_program,\"lintwopi\");\n",
       "        gl.enableVertexAttribArray(ltp);\n",
       "        gl.vertexAttribPointer(ltp, 1, gl.FLOAT, 0, 0,0); // 4 = size of float\n",
       "        gl.uniformMatrix4fv(reboundView.orbit_shader_mvp_location,false,mattransp(mvp));\n",
       "\n",
       "        // Need to do this one by one\n",
       "        // because WebGL is not supporting\n",
       "        // instancing:\n",
       "        for(i=0;i<reboundView.N-1;i++){\n",
       "            var focus = new Float32Array(reboundView.orbit_data.buffer,4*9*i,3);\n",
       "            gl.uniform3fv(reboundView.orbit_shader_focus_location,focus);\n",
       "            var aef = new Float32Array(reboundView.orbit_data.buffer,4*(9*i+3),3);\n",
       "            gl.uniform3fv(reboundView.orbit_shader_aef_location,aef);\n",
       "            var omegaOmegainc = new Float32Array(reboundView.orbit_data.buffer,4*(9*i+6),3);\n",
       "            gl.uniform3fv(reboundView.orbit_shader_omegaOmegainc_location,omegaOmegainc);\n",
       "\n",
       "            gl.drawArrays(gl.LINE_STRIP,0,500);\n",
       "        }\n",
       "    }\n",
       "}\n",
       "require.undef('rebound');\n",
       "define('rebound', [\"jupyter-js-widgets\"], function(widgets) {\n",
       "    var ReboundView = widgets.DOMWidgetView.extend({\n",
       "        render: function() {\n",
       "            this.el.innerHTML = '<canvas style=\"display: inline\" id=\"reboundcanvas-'+this.id+'\" style=\"border: none;\" width=\"'+this.model.get(\"width\")+'\" height=\"'+this.model.get(\"height\")+'\"></canvas>';\n",
       "            this.model.on('change:t', this.trigger_refresh, this);\n",
       "            this.model.on('change:count', this.trigger_refresh, this);\n",
       "            this.startCount = 0;\n",
       "            this.gl = null;\n",
       "            // Only copy those once\n",
       "            this.scale = this.model.get(\"scale\");\n",
       "            this.width = this.model.get(\"width\");\n",
       "            this.height = this.model.get(\"height\");\n",
       "            this.orbits = this.model.get(\"orbits\");\n",
       "            this.orientation = this.model.get(\"orientation\");\n",
       "            startGL(this);\n",
       "        },\n",
       "        trigger_refresh: function() {\n",
       "            updateRenderData(this);\n",
       "            drawGL(this);\n",
       "        },\n",
       "    });\n",
       "    return {\n",
       "        ReboundView: ReboundView\n",
       "    };\n",
       "});\n",
       "      \n",
       "</script>\n"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "67b74bd5df21427abf8dd3777be3c608"
      }
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sim.getWidget(size=(400,400))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "import time\n",
    "\n",
    "time.sleep(10)\n",
    "\n",
    "for t in np.arange(0., 40., 0.1):\n",
    "    sim.integrate(t)\n",
    "    time.sleep(0.1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "Finally, I did a quicktime screencast and then converted that into a gif via the following command:\n",
    "\n",
    "ffmpeg -i Leela\\'s\\ Constellation\\ animation.mov -pix_fmt rgb24 -r 10 -f gif - | gifsicle --optimize=3 --delay=3 > Leela\\'s\\ Constellation\\ animation.gif"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}