Skip to content

Instantly share code, notes, and snippets.

@AlexanderFabisch
Created May 2, 2020 13:34
Show Gist options
  • Save AlexanderFabisch/82e49799705e223906df0311cc4dfbe4 to your computer and use it in GitHub Desktop.
Save AlexanderFabisch/82e49799705e223906df0311cc4dfbe4 to your computer and use it in GitHub Desktop.
pytransform3d issue 43
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Populating the interactive namespace from numpy and matplotlib\n"
]
}
],
"source": [
"%pylab inline\n",
"from pytransform3d.rotations import *"
]
},
{
"cell_type": "code",
"execution_count": 25,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(1.31907686094078e-17, 0.0, 4.440892098500626e-16)"
]
},
"execution_count": 25,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAsgAAAFgCAYAAACmDI9oAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdeXxW5Z3//9cnCSGETQhh35VFVlG2alVsq1Jr1SqOWKe2rmOtdvlOZ+osXaYdx87Mr9OOW61tqcu0LqW2tdXWqh3EBSuguCCriBKQfQtLCEmu3x+JzC0kECB3Esjr+XjE5Fz3Odf53DmPkLdXrnOdSCkhSZIkqVpOUxcgSZIkNScGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAx5TV1AY+jSpUvq379/U5chSZKkZmTu3LnrU0rFe7e3iIDcv39/5syZ09RlSJIkqRmJiHdqa3eKhSRJkpTBgCxJkiRlMCBLkiRJGVrEHGRJkqQjwe7duykpKaGsrKypSzmqFBQU0Lt3b1q1alWv/Q3IkiRJzURJSQnt27enf//+RERTl3NUSCmxYcMGSkpKGDBgQL2OcYqFJElSM1FWVkZRUZHhuAFFBEVFRQc1Km9AliRJakYMxw3vYL+nR+QUi4hoC9wJlAMzUko/b+KSJEmSdJRoNiPIETEtItZGxBt7tU+OiEURsTQibqppvhCYnlK6Bjiv0YuVJElqQVatWsWUKVOauoxG02wCMnAPMDmzISJygTuAjwPDgEsjYhjQG1hRs1tlI9YoSZLU4vTs2ZPp06c3eL+Vlc0zxjWbgJxSmgls3Kt5PLA0pbQspVQOPAicD5RQHZKhGb0HSZKkI9nXvvY17rzzzj3b3/rWt/je977H8uXLGTFiBADLly/n1FNP5cQTT+TEE0/khRdeqLWvCy64gJNOOonhw4dz991372lv164d3/jGN5gwYQKzZs1i9uzZnHzyyYwePZrx48dTWlpKZWUlX/3qVxk5ciSjRo3itttuA+Dpp59mzJgxjBw5kiuvvJJdu3YBcNNNNzFs2DBGjRrFV7/61cP+PjT3Oci9+L+RYqgOxhOAW4HbI+ITwO9qOzAirgWuBejbt2+Wy5QkSWpY//K7+by5amuD9jmsZwe++cnhdb4+depUvvzlL3P99dcD8PDDD/PHP/6RqqqqPft07dqVJ598koKCApYsWcKll17KnDlz9ulr2rRpdO7cmZ07dzJu3DguuugiioqK2L59OyNGjODb3/425eXlDB06lIceeohx48axdetW2rRpw913383bb7/NK6+8Ql5eHhs3bqSsrIzPfe5zPP300wwePJjLL7+cH/7wh1x++eX8+te/ZuHChUQEmzdvPuzvU3Mffa3tlsOUUtqeUroipfT5um7QSyndnVIam1IaW1xcnOUyJUmSjnxjxoxh7dq1rFq1ildffZVOnTrtM9C4e/durrnmGkaOHMnFF1/Mm2++WWtft956K6NHj2bixImsWLGCJUuWAJCbm8tFF10EwKJFi+jRowfjxo0DoEOHDuTl5fHUU09x3XXXkZdXPZbbuXNnFi1axIABAxg8eDAAn/3sZ5k5cyYdOnSgoKCAq6++mkceeYTCwsLD/j409xHkEqBPxnZvYFUT1SJJktRo9jfSm01Tpkxh+vTprF69mqlTp+7z+ve//326devGq6++SlVVFQUFBfvsM2PGDJ566ilmzZpFYWEhkyZN2rMOcUFBAbm5uUD1QzxqW4KttvaUUq315uXl8dJLL/H000/z4IMPcvvtt/PnP//5oN93puY+gjwbGBQRAyIiH5gKPNrENUmSJB21pk6dyoMPPsj06dNrXbliy5Yt9OjRg5ycHO6///5ab7TbsmULnTp1orCwkIULF/Liiy/Weq6hQ4eyatUqZs+eDUBpaSkVFRWcddZZ3HXXXVRUVACwceNGhg4dyvLly1m6dCkA999/P6effjrbtm1jy5YtnHPOOfzgBz9g3rx5h/09aDYBOSIeAGYBQyKiJCKuSilVADcATwALgIdTSvObsk5JkqSj2fDhwyktLaVXr1706NFjn9evv/567r33XiZOnMjixYtp27btPvtMnjyZiooKRo0axde//nUmTpxY67ny8/N56KGHuPHGGxk9ejRnnnkmZWVlXH311fTt25dRo0YxevRofvGLX1BQUMDPfvYzLr74YkaOHElOTg7XXXcdpaWlnHvuuYwaNYrTTz+d73//+4f9PYi6hquPJmPHjk21TR6XJElqThYsWMDxxx/f1GUclWr73kbE3JTS2L33bTYjyJIkSVJzYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmStF+rVq2q9aEhje3kk09ulPMYkCVJkrRfPXv2ZPr06Q3eb21P4dufF154ocFrqI0BWZIkSQB87Wtf484779yz/a1vfYvvfe97LF++nBEjRgCwfPlyTj31VE488UROPPHEOkPrBRdcwEknncTw4cO5++6797S3a9eOb3zjG0yYMIFZs2Yxe/ZsTj75ZEaPHs348eMpLS1l/vz5jB8/nhNOOIFRo0axZMmSPccCzJgxg0mTJjFlyhSGDh3KZZddRkM+/C6vwXqSJElSw/nDTbD69Ybts/tI+Ph363x56tSpfPnLX+b6668H4OGHH+aPf/wjVVVVe/bp2rUrTz75JAUFBSxZsoRLL72U2p5YPG3aNDp37szOnTsZN24cF110EUVFRWzfvp0RI0bw7W9/m/LycoYOHcpDDz3EuHHj2Lp1K23atOGuu+7iS1/6Epdddhnl5eW1jjS/8sorzJ8/n549e3LKKafw/PPP8+EPf7gBvkkGZEmSJNUYM2YMa9euZdWqVaxbt45OnTrRt29fli9fvmef3bt3c8MNNzBv3jxyc3NZvHhxrX3deuut/PrXvwZgxYoVLFmyhKKiInJzc7nooosAWLRoET169GDcuHEAdOjQAYAPfehD3HzzzZSUlHDhhRcyaNCgffofP348vXv3BuCEE05g+fLlBmRJkqSj2n5GerNpypQpTJ8+ndWrVzN16tR9Xv/+979Pt27dePXVV6mqqqKgoGCffWbMmMFTTz3FrFmzKCwsZNKkSZSVlQFQUFBAbm4uACklImKf4z/96U8zYcIEHnvsMc4++2x+8pOf8JGPfOQD+7Ru3XrP17m5uVRUVBzW+87kHGRJkiTtMXXqVB588EGmT59e68oVW7ZsoUePHuTk5HD//ffXOv1hy5YtdOrUicLCQhYuXMiLL75Y67mGDh3KqlWrmD17NgClpaVUVFSwbNkyBg4cyBe/+EXOO+88XnvttYZ9kwdgQJYkSdIew4cPp7S0lF69etGjR499Xr/++uu59957mThxIosXL6Zt27b77DN58mQqKioYNWoUX//615k4cWKt58rPz+ehhx7ixhtvZPTo0Zx55pmUlZXx0EMPMWLECE444QQWLlzI5Zdf3uDvc3+iIe/4a67Gjh2baps8LkmS1JwsWLCA448/vqnLOCrV9r2NiLkppbF77+sIsiRJkpTBgCxJkiRlMCBLkiQ1Iy1h+mtjO9jvqQFZkiSpmSgoKGDDhg2G5AaUUmLDhg21LkdXF9dBliRJaiZ69+5NSUkJ69ata+pSjioFBQV7HipSHwZkSZKkZqJVq1YMGDCgqcto8ZxiIUmSJGUwIEuSJEkZDMiSJElSBgOyJEmSlMGALEmSJGUwIEuSJEkZjsiAHBHHR8RdETE9Ij7f1PVIkiTp6NHoATkipkXE2oh4Y6/2yRGxKCKWRsRN++sjpbQgpXQd8FfA2GzWK0mSpJalKUaQ7wEmZzZERC5wB/BxYBhwaUQMi4iREfH7vT661hxzHvAc8HTjli9JkqSjWaM/SS+lNDMi+u/VPB5YmlJaBhARDwLnp5RuAc6to59HgUcj4jHgF3u/HhHXAtcC9O3bt8HqlyRJ0tGtuTxquhewImO7BJhQ184RMQm4EGgNPF7bPimlu4G7AcaOHZsaqlBJkiQd3ZpLQI5a2uoMtSmlGcCMbBUjSZKklqu5rGJRAvTJ2O4NrGqiWiRJktSCNZeAPBsYFBEDIiIfmAo82sQ1SZIkqQVqimXeHgBmAUMioiQirkopVQA3AE8AC4CHU0rzG7s2SZIkqSlWsbi0jvbHqeOGO0mSJKmxNJcpFpIkSVKzYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMuQ1dQGHIiJygO8AHYA5KaV7m7gkSZIkHSUafQQ5IqZFxNqIeGOv9skRsSgilkbETQfo5nygF7AbKMlWrZIkSWp5mmIE+R7gduC+9xsiIhe4AziT6sA7OyIeBXKBW/Y6/kpgCDArpfSjiJgOPN0IdUuSJKkFaPSAnFKaGRH992oeDyxNKS0DiIgHgfNTSrcA5+7dR0SUAOU1m5W1nScirgWuBejbt2+D1C5JkqSjX3O5Sa8XsCJju6SmrS6PAGdHxG3AzNp2SCndnVIam1IaW1xc3HCVSpIk6ajWXG7Si1raUl07p5R2AFdlrxxJkiS1VM1lBLkE6JOx3RtY1US1SJIkqQVrLgF5NjAoIgZERD4wFXi0iWuSJElSC9QUy7w9AMwChkRESURclVKqAG4AngAWAA+nlOY3dm2SJElSU6xicWkd7Y8DjzdyOZIkSdIHNJcpFpIkSVKzYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKkNfUBRytxv7rk2zasbupy2g00VD91NFR1HWGejRn9vl+Px9se3+/+OCx8X+f9rwWe23v6Sv2vJYT8YH9Iqr3yYn4wOs5ER/8Ouf914LcgNyc2LNPXm7N55wgN+eD23m5ObTKDfJycsjLDVrVbLfKzSE/L4f8jM+tW+VQkJdLQatcWrfKoU2r6q/b5udR2DqXdq3zaJ2Xs+f9ZUtVVaK8soqCVrlZPY8kSYfCgJwlV5wygJ3llU1dRqNIpIbpp45u6uq97v1T5sY+X6aMA9//Mu2znfY5R0rVrZn7pERGW812gqqMfffslxJVNftXVSUSiaqq6n2ral6vrPm6qipRWZWoSomKqirKKlJ1W0pUVFa/VlmVqKhKVFRWsbvmc0VlYndVFbtr9jkUuTlBYX51aO7YphUdC1vRsU0rjmlT87mwFV3bF1DcoTXd2hfQtUNrOhfmk5Oz/1CdUmLBe6X8dt5KHn11Fau3ljG4a3tO6HMMY/oewwl9j2FQ1/bkHqAfSZKyLVJdKeMoMnbs2DRnzpymLkNqVJVVid2VVZRXVlFeUf2xq6KKst2VNR9VlFVUUlZeyY7ySnaUV7BtVyXbd1WwvbyCbWUVbC3bzeYdu9mys/pj847d7Ny97//45eUE3ToU0L9LIQO7tGNAl7YMLG7LscXtAHj01VX8dt5KFq/ZRl5OcPrgYob17MDrK7cwb8VmNtf8taVtfi4fOraIS8b15YwhxeTlOgtMkpQ9ETE3pTR273ZHkKWjVG5OkJuT2+DTGMp2V7KudBdrS8tYs3UXa7eWsaZ0F+9t3snb67fzm1dWUrqrYp/jxvbrxHcuGMEnRvagc9v8Pe0pJZZv2MEr727ilXc388T81Ty1YA7dOxRwybg+XDKuDz2PadOg70GSpP1xBFlSg0opsX5bOcvWbePt9dvZXl7JWcO60adzYb2O311ZxdML1vKLl97l2SXrCOCMIV2ZMLDzB6a0vD91ZViPDpw6qIujzZKkg1bXCLIBWVKztWLjDh6c/S4PzylhXemuOvfr1qE1F53Ym4vH9mFAl7aNWKEk6UhmQDYgS0esyqrEzt2Ve1YNiZpVQ6pSYubidTw8p4QZi9ZSlWD8gM781dg+nDOyO4X5ziKTJNXNgGxAlo5qa7aW8auXS/jlnBLeXr+d9q3zOO+Enlw6vi8jenVs6vIkSc3QQQfkiLiwHv2WpZQeP9ziss2ALLUcKSVmL9/Egy+9y2Ovv8euiipG9OrA1HF9Of+EnrQvaNXUJUqSmolDCcgbgN+y/2dAnJZSOrZhSsweA7LUMm3ZsZvfzFvJAy+9y8LVpRS0yqFHxza0rXkoSrvWebRtnUdhfu6eJfB21ix7V7a7ksL8PP79olH0LarfDYaSpCPLoQTk/0kp/fUBOj3gPs2BAVlq2VJKvFqyhUfnrWJtaVn1Ws+7KtlWs+bz9l2VtM7LoU1+LoX5ubRplUub/FzmvrOJfkWFTL/uZJ/6J0lHIecgG5AlHaSn3lzD1ffN4dMT+vJvnxq5331TSixbv53l67ezfMMO3tlQ/Xn5+u307VzIHZ8+kY6FTu+QpObkoB8UcqA5yCmlRxqiMElqrj42rBufn3QsP5zxFuP6d+JTY3rXul9p2W5ufOAVZixat6etfUEeA7q0ZUSvDjz55ho++7OX+J+rJ9CutStrSFJzt79/qT9Z87krcDLw55rtM4AZgAFZ0lHvb88czMvvbOIfH3mD4T07Mrhb+w+8vmLjDq66dzbL1m3n7ycPYeLAIvoXtaVTYSsiqm/heGL+aq7/+ctcec9s7r1iPG3yna4hSc1ZnY+eSildkVK6guoHVw1LKV2UUroIGN5o1dUhIoZFxMMR8cOImNLU9Ug6euXl5nDbpWNo2zqP6/5nLtsyHqM9951NXHDH86zeUsa9V47n+knHcWLfTnRum78nHAOcPbw737/kBOYs38i198+hbHdlU7wVSVI91efZrP1TSu9lbK8BBh/qCSNiWkSsjYg39mqfHBGLImJpRNx0gG4+DtyWUvo8cPmh1iJJ9dG1QwG3XTqG5eu38w+PvE5Kid/OW8mlP36RdgV5PHL9KZxyXJf99nHe6J78+0WjeHbJer7w85cpr6hqpOolSQerPpPhZkTEE8ADVI8mTwX+9zDOeQ9wO3Df+w0RkQvcAZwJlACzI+JRIBe4Za/jrwTuB74ZEecBRYdRiyTVy4eOLeKrZw/hP/64iC07dzNz8TrGD+jMj/76JDq1za9XHxeP7UNZRRVf/80bfPmhV7h16hjycuszTiFJakwHDMgppRtqbtg7tabp7pTSrw/1hCmlmRHRf6/m8cDSlNIygIh4EDg/pXQLcG4dXX2hJljXOhc6Iq4FrgXo27fvoZYrSXtcd9qxzF2+iacXrmXKSb35t0+NJD/v4ALuZyb2Y9fuSv71sQUsXfss/YraUty+NV3ataa4fWuK2+Uzus8x9OjYJkvvQpJ0IPW6nbpmxYps3pTXC1iRsV0CTKhr55qA/Y9AW+A/a9snpXQ3cDdUL/PWQHVKasFycoLbPj2Gee9u5kPHFn1gnvHBuPrUgbTJz+UPr69mxcYdvPzOJjbuKOf9VTdzAk4fXMzU8X35yNCutHKUWZIa1QEDckRMBG4DjgfyqZ72sD2l1KEB66jtt0ydoTaltJya0WFJakyF+XmcfID5xvVx2YR+XDah357tisoqNm4vZ83WXfzpzdU8PGcFf3P/XIrbt+bik3pzybg+9Ctqe9jnlSQdWH1GkG+net7xL4GxVN8Ud1wD11EC9MnY7g2sauBzSFKzlZebQ9cOBXTtUMDI3h350kcH8b+L1vHgS+9y1zNvceeMt7h+0rH8/eShTV2qJB316jvFYmlE5KaUKoGfRcQLDVzHbGBQRAwAVlIdyD/dwOeQpCNGXm4OZw7rxpnDuvHelp1870+LuXPGW3Rs04q/Of3Y/R5btruS/3nxHdaW7mL7rgp2lFeyo7z6c25O8LXJQzm+R0P+EVCSji71Ccg7IiIfmBcR/wG8R/Xc30MSEQ8Ak4AuEVECfDOl9NOIuAF4guopHNNSSvMP9RySdDTp0bEN/3HRKMp2V3LLHxZS1K41U06q/al+60p3cc19c5i3YjNtWuXStnUubfJzaZufR2F+Lu9u3MFf/WgWP758LBMHugiQJNUmUtr//WsR0Y/qtY/zga8AHYE7U0pLs19ewxg7dmyaM2dOU5chSYdlV0UlV90zh1nLNnD3Z07io8d3+8Dri9eUcsXPZrNh+y5+cMkYJo/ovk8fqzbv5PJpL/Huhh38YOoJnDOyR2OVL0nNTkTMTSmN3bt9v7dG1yyjdnNKqSyltDWl9C8ppf93JIVjSTpatM7L5a7PnMTwnh34wi9eZu47G/e89szidVx05wvsrqzi4b/5UK3hGKDnMW2Yft2HGNm7I1/4xcvcN2t5rfu9vX473/3DQr7y0Dx2VfjkP0kty36nWKSUKiOiOCLyU0rljVWUJKl27Vrn8bPPjePiu2Zxxc9m88vrTmb28o1889H5DOrajmmfG0fPY/a/hvIxhfn8/OoJ3PCLV/jGb+ezZmsZXz1rCLsqqvjDG+/x4Esr+MvbG8nNCSqrEv2L2vKljw1qpHcoSU2vPlMsfgScCDwKbH+/PaX0X9ktreE4xULS0WbFxh1MuesFtpVVsL28ko8M7cqtl46hXet63XsNVC8t9/XfzueBl95lfP/OLFi9ldKyCvoVFXLJuD5MObE3//rYAv7wxns8/sVTGdSt/X77e3fDDjbvLGdU72MO9+1JUqOoa4pFff4lXVXzkQPs/19HSVKj6NO5kHuvHM8VP5vNJeP68k+fOJ7cnIN7cElebg7/9qkRdOvQmh89s4yzh3fjknF9mTCgMzk1fX3jk8OYuWQdX/vVa0y/7uQ97Xt7d8MOPnXn82zYXs6pg7rw5Y8N5qR+nQ77fUpSUzjgCPLRwBFkSUerlNIhP9Gvvv088nIJ/+/hV/n2+cO5/EP993l9845yLvzhC2zYVs7nTu7P/S++w8bt5UwaUsxXPjaY0X0cUZbUPB30TXoR8a16dHrAfSRJ2dMQ4fhA/XxqTC9OHdSFf//DQlZu3vmB13ZVVHLt/XMp2biTH18+lq+cOZhn//4MvjZ5KPNWbOb8O57n6ntn8+KyDVRUVjVIrZKUbXWOINesUby/ecYBXJNSavaPdXIEWZIOz4qNOzjr+zOZOLAz0z43joigqirx5Yfm8eirq7j10jGcN7rnB44pLdvNvS8s5+6Zy9haVkGnwlZ89Pjqh5+cNqiYNvm5TfRuJKnaocxB/jEHnnP848OqSpJ0ROjTuZCvnj2E7/z+TX732nucN7on33tyEY++uoq/nzxkn3AM0L6gFTd8ZBBXnDKAGYvW8eSbq3li/mqmzy2hoFUOHz6umPNO6MnZw7vROs+wLKn5cA6yJKleKqsSF975PCWbdnLNaQP57h8Wcun4Pvzbp0bWe6rH7soqXnp7I3+av5o/vbmG97aU0amwFVNO6s3U8X05trhdlt+FJP2fukaQDciSpHpb8N5WPnnbc1RUJU4fXMxPPzuWvNz9PnOqTlVVieeWrueBl97lyTfXUFGVmDCgM5eO70u/okKqUqKyqjqYV9X8rhrWowOd2uY35FuS1IIZkA3IktQgfvrc28xcvI47LjvxoNZd3p+1pWVMn1vCgy+t4N2NO+rcLwJG9z6G0wYXc/rgYkb37njIAV2SDikg1zxq+osppe9ns7hsMyBL0pGhqirxyopNbN1ZQU5OkBtBTg7k5eSwu7KK2cs38szidby6YjNVCToU5PHhQV04sW8nRvc5huE9O1CY3zChXdLR75BHkCNiRkppUrYKawwGZEk6umzeUc7zSzfwzOK1PL90w57l53ICBndrz+jexzCid0d6dCiga4fWFLdvTZd2rWnlaLOkDIcTkG8GOgIP8cFHTb/c0EVmiwFZko5u60p38VrJZl5dsZl5JVt4rWQzm3fs3me/zm3zGdilLX939hAmDCxqgkolNSeHE5D/t5bmlFL6SEMVl20GZElqWVJKvLeljLWlu1hX87G2tIx1pbuYsWgdKzfv5MIxvfiHc46nuH3rOvtZV7qLqpTo1qGgEauX1FgOZR1kAFJKZ2SnJEmSsiMi6HlMG3oe02af13aWV3L7/y7h7pnLeHLBGr561hD+emI/cnOql6pbV7qLP77xHo+9/h5/eXsjORFMObE3N370OHp3KmzstyKpCdRnBLkj8E3gtJqmZ4Bvp5S2ZLm2BuMIsiRpb2+t28Y3fzuf55auZ0SvDlxwQi+eXrCWv7y9gaoExxa35dxRPdmycze/+Mu7JBKXju/LF844zhFl6ShxOFMsfgW8Adxb0/QZYHRK6cIGrzJLDMiSpNqklPj9a+/xnd+/ydrSXQzs0pZzR/XgE6N6Mrhbuz0PQFm1eSe3/Xkpv5yzgtyc4DMT+3Ht6QPp2t6gLB3JDicgz0spnXCgtubMgCxJ2p8d5RWs3bqLfkWF+30q4LsbdvDfTy/h16+UAHDKcV345KienD28Ox0LWzVWuZIayOEE5FnA36WUnqvZPgX4/1JKH8pKpVlgQJYkNaRl67bxyMsr+d1rq3hnww5a5QanD+7KJ0f34Kxh3WmTn9vUJUqqh8MJyKOB+6he6g1gE/DZlNJrDV5llhiQJUnZkFLi9ZVb+N2rq/j9a+/x3pYyenYs4OvnDmPyiO77HY2W1PQO9Ul6OcCUlNLDEdEBIKW0NXtlZocBWZKUbVVViRfe2sDNjy9gwXtbOXVQF/7lvOEMLG5X5zE7yyspaJVjkJaayOGMIM9MKZ22352aOQOyJKmxVFRW8T8vvsP3/rSYXRVVXHPaAG44YxBt8nNZv20XL729kb8s28CLyzayaE0p3TsUcPKxRZx8XBdOOa6IHh33XZpOUnYcTkD+OrCTfZ+kt7Ghi8wWA7IkqbGtK93FLX9YwCMvr6RnxwLats5jydptABTm53JSv06M6XMMy9Zv54W3NrBxezkAA7u05eTjijhtUDEnH9eFdq0P+MgCSYfocALy27U0p5TSwIYqLtsMyJKkpvLS2xv5rycXUdAqlwkDipg4sDMjenWkVW7Onn2qqs37YM8AAB1xSURBVBILV5fywlvreX7pev7y9kZ2lFeSlxOc1K8Tpw8p5vTBxQzr0cHpGFIDOpw5yB9KKT2fzeKyzYAsSTqSlFdUMeedjcxcvJ5nFq9jwXvVt/90adeacf07cVK/6o/hPTuSn5dzgN4k1eWwlnlr6iXdImIg8E9Ax5TSlJq2C4BPAF2BO1JKf6rreAOyJOlItnZrGTOXrOe5JeuY++4mVmzcCUDrvBxG9zmGU47twrWnDXR5OekgHU5A/hfgNeCRdKCdaz9+GnAusDalNCKjfTLw30Au8JOU0nfr0df09wNyRlsnqtdlvqqu4wzIkqSjyZqtZcx9ZxNz39nEnHc28VrJZo4tbsdtl47h+B4dmro86YhxOAG5FGgLVFJ9s15QPQe5Xj+BEXEasA247/2AHBG5wGLgTKAEmA1cSnVYvmWvLq5MKa2tOa62gPw94OcppZfrqsGALEk6mj2/dD1feWgem3fu5p8/cTyfmdjPucpSPdQVkA84cSml1D6llJNSapVS6lCzXe//PU0pzQT2XvFiPLA0pbQspVQOPAicn1J6PaV07l4fa+t4QxER/w78obZwHBHXRsSciJizbt26+pYrSdIR55TjuvCHL53KKccW8Y3fzuea++ayqWZVDEkH74ABuSaI/nXNcm9ERJ+IGH+Y5+0FrMjYLqlpq6uGooi4CxgTEf9Q03wj8DFgSkRct/cxKaW7U0pjU0pji4uLD7NcSZKat6J2rZn2uXF8/dxhPLN4LR//72d5bsl6DmF2pNTi1WdxxTuBKuAjwHeoni5xBzDuMM5b29996vwJTiltAK7bq+1W4NbDqEGSpKNKRHDVhwcwYUBnbnzgFf76p39hSLf2XDy2NxeM6UWXdq2bukTpiFCftWEmpJS+AJQBpJQ2AfmHed4SoE/Gdm9g1WH2KUmSgBG9OvLYFz/MzZ8aQZv8XP71sQVM/Lenuea+Ofxp/mp2V1YddJ9VVYmy3ZVs21XBpu3lrCvd5ei0jlr1GUHeXXNTXQKIiGKqR5QPx2xgUEQMAFYCU4FPH2afkiSpRmF+HpdN6MdlE/qxZE0p0+eW8KuXV/Lkm2to0yqXwd3aMbhbe4Z0b8/Q7h0Y0r09uTnB4jWlLFlTyuI121i8ppSla7exeeduKqv2DcOjenfk3z41khG9OjbBO5Sypz6rWFwGXAKcCNwLTAH+OaX0y3qdIOIBYBLQBVgDfDOl9NOIOAf4AdUrV0xLKd18qG/iQFzFQpIk2F1ZxTOL1vHCWxtYtGYri1ZvY/22XbXu2751HoO6tWNQ1/YUt29Nq9wcWuUFrXJyaJUb7Kqo4ifPvc2Gbbv47Mn9+duzhvhYbB1xDnmZt5qDhwIfpXru8NMppQUNX2L2GJAlSardhm27WLS6lIWrS6lKiUHd2jO4Wzu6dyg44FJxW3bu5j+fWMjP//Iu3doX8K3zhnH28O4uMacjxmEF5COdAVmSpOx5+d1N/OMjr7NwdSkfHdqVmz4+lEHd2jd1WdIBGZANyJIkZc3uyip+9vzbfP/JJezcXclpg4u54pT+nD6omJwcR5TVPBmQDciSJGXdhm27eOCld7lv1jusLd3FwOK2XHHKAC46sReF+c5RVvNyuHOQ+wGDUkpPRUQbIC+lVJqFOrPCgCxJUuMqr6ji8dffY9rzb/NayRbat87j+B4d6N25Db07FdKnUxv6dC5kYJe2dO1Q0NTlqoWqKyAf8H/lIuIa4FqgM3As1WsW30X1TXuSJEn7yM/L4YIxvTj/hJ68/O4mfjmnhGXrtvPiWxt4b+tKMsfnvviR4/jyxwY7FUPNRn3+1vEFYDzwF4CU0pKI6JrVqiRJ0lEhIjipX2dO6td5T1t5RRWrNu9kxaYd/PqVldz656UsXF3Kf11yQp1LxaWUeGrBWl4r2czFJ/Whb1FhY70FtUD1Cci7Ukrl7y/ZEhF57Oex0JIkSfuTn5dD/y5t6d+lLR8+rgsjenbk5scXcOGdz/Pjy8fSr6jtB/Z/Y+UWbn5sAbOWbQDgzhlvccEJvfjCGccysLhdU7wFHeXq86jpZyLiH4E2EXEm8Evgd9ktS5IktQQRwZUfHsC9V4xnzdZdnH/H8zy/dD0Aa7aW8Xe/fJVP3v4cC1dv5TvnD+f5mz7CZz/Un8deX8XH/usZvvjAKyxec8TcFqUjRH2epJcDXAWcRfWDQp4AfpKOoOUvvElPkqTm750N27nmvjm8tW47543uyR/fWE1FVRVXnDKAL5xxHB3btNqz77rSXfzkuWXcP+sddu6uZFy/zhS3b02HNq3o2KYVHdrk0bFNK4Z278CJfY/x4SWqlcu8GZAlSWr2tu2q4CsPzePJN9dwzsjufG3y0H2mXGTauL2cac+9zfNvrWfLzt1s3bmbLTt3s7vy//LNyF4dufLD/fnEyJ7k59Xnj+dqKQ46IEfE6+xnrnFKaVTDlZddBmRJko4cKSVWby2jR8c2h3x82e4qNu8s538XrmPa82+zdO02itu35vKJ/fj0hL4UtWvdwFXrSHQoAbnf/jpMKb3TQLVlnQFZkqSWK6XEs0vWM+35t5mxaB35eTl8fER3zhnZg9MHF1PQKrepS1QTOeh1kI+kACxJklSXiOC0wcWcNriYpWu3ce8Ly/n9a6v47bxVFObn8pGhXfnEyB5MGtKVNvkHDss7yitYvGYb727cwYeP60LntvmN8C7UmOpzk14p+0612ALMAf42pbQsS7U1GEeQJUlSpt2VVfxl2UYef+M9nnhjNRu2l9OmVS6DurWjc9t8itq2pqhdPp3b5tOpsBWrNpexcPVWFq0u5Z2NO/Y86KRDQR5/e9YQLpvQl7xc5zcfaQ75Jr2I+BdgFfALqlexmAp0BxYBn08pTWrwahuYAVmSJNWlorKKl5Zv5Ik3VrN8ww42bi9nw7ZdrN9eTnlFFQA5Af2L2jKke3uGdu/AkO7t6dw2n/9+ejHPL93AkG7t+eYnh3HycV2a+N3oYBxOQP5LSmnCXm0vppQmRsSrKaXRDVxrgzMgS5Kkg5VSYnt5JZu2l9OlXetap1+klHhi/hr+9bE3Kdm0k3NGducfzzme3p180t+R4KDnIGeoioi/AqbXbE/JeO3oXyNOkiS1SBFBu9Z5dT7++v19Jo/ozqQhxfx45jLumLGUpxas5ZRjizh9cDGnD+lK/6JC12E+wtRnBHkg8N/Ah6gOxC8CXwFWAiellJ7LdpGHyxFkSZLUGFZu3smPZy7jmcXreHv9dgD6di7k9MHFnDG0mNMHdyU3x7DcXPigEAOyJElqRO9s2M7Mxet4ZvE6XnhrAzvKKxlY3JYbzjiO80b39Ka+ZuBw5iAXA9cA/cmYkpFSurKBa8waA7IkSWpKuyoqefLNNdz+56UsXF1Kv6JCrp90LJ8a09un+zWhwwnILwDPAnOByvfbU0q/augis8WALEmSmoOqqsSTC9Zw25+X8MbKrfQ6pg1XnNKfk/p1Ykj39hTm1+f2MDWUwwnI81JKJ2StskZgQJYkSc1JSokZi9Zx65+X8Mq7mwGIgH6dCxnavQNDe7RneM+OjOl7DF18LHbWHM4qFr+PiHNSSo9noS5JkqQWJyI4Y2hXJg0pZsXGnSxYvZWF75WycPVWFq4u5Yk3V+95GEm/okJO6tuJE/t14sS+nRhY3JbaFsXIz81xtYwGUt8n6bUFdgG7qX5YSEopdch+eQ3DEWRJknQk2VFewfxVW3n5nU28/O4m5r6zmfXbdu33mFG9O3LbpWPoV9S2kao88rmKhQFZkiQdoVJKlGzaydx3NrFy8859Xi+vqOKeF5ZTVZX4z4tHMXlEjyao8shzOFMsiIhOwCCg4P22lNLMhitPkiRJdYkI+nQupE/nup/QN+Wk3tzwi5e57n9e5spTBnDTx4e6QsYhOuB3LSKuBmYCTwD/UvP5W9kt6wPnHxgRP42I6RltkyLi2Yi4KyImNVYtkiRJzVWfzoX88rqT+dzJ/Zn2/Nv81Y9m1TrarAOrzwjyl4BxwIsppTMiYijVQfmAImIacC6wNqU0IqN9MtVP58sFfpJS+m5dfaSUlgFXZQZkqp/ot43qEe2S+tQiSZJ0tMvPy+Fb5w1n/IDO/P301/jErc9ywQm9yNnr5r2I6jnLZw3rTpv83Dr7Synx8rubeHbJek7q14mTj+3SIp4EWJ+AXJZSKosIIqJ1SmlhRAypZ//3ALcD973fEBG5wB3AmVSH29kR8SjVYfmWvY6/MqW0tpZ+n00pPRMR3YD/Ai6rZz2SJElHvXNG9uD4Hh3424fn8au5+44l7q6q4qfPVdGudR6TR3TnwjG9mDiwiJya8Lt8/XZ+/cpKfjNvJe9s2LHnuG4dWnP+Cb341JheHN/jiFmv4aDVJyCXRMQxwG+AJyNiE7CqPp2nlGZGRP+9mscDS2tGhomIB4HzU0q3UD3aXJ9+q2q+3AS4OKAkSdJeBnRpyyPXn1Lra1VViRff3sCvX17JH95YzfS5JfToWMBZw7rx+sotvPzuZiLg5GOLuPEjgzhjSDF/eXsjj7y8kmnPvc3dM5cxtHt7LjqxN389sd9+R6GPRAe1ikVEnA50BP6YUiqv5zH9gd+/P8UiIqYAk1NKV9dsfwaYkFK6oY7ji4CbqR5x/klK6ZaIuBA4GzgG+GFKaUYtx10LXAvQt2/fk9555516v09JkqSWomx39WOwf/3KSp5ZvI7jitvxqRN7cf4JPenRsc0++2/cXs7vX1vFIy+vZN6KzfTu1IbvnD+CM4Z2bYLqD0+TLfNWS0C+GDh7r4A8PqV0Y7ZqcJk3SZKkAyuvqKJVbtT7gSMvLtvAP/36dd5at51zRnbnm58cTrcOBQc+sJmoKyA3xdofJUCfjO3e1HPKhiRJkrInP+/gnsY3cWARj3/pVL561mCeXrCWj37vGe55/m0qq47s52zUax3kBjYbGBQRA4CVwFTg001QhyRJkg6kfAc8eiOcfTO0777Py63zcrnhI4P45Oie/PNv3uBbv3uTe15YTpd2H7xN7PjyN7ii/OcM7FzLCHPbLjD159l6BwctqyPIEfEAMAsYEhElEXFVSqkCuIHq9ZQXAA+nlOZnsw5JkiQdovmPwBvT4elv73e3fkVtue/K8dx66Rj6FrWldaucD3wcW/oSA7fPg7zWtXw0r2kZWR1BTildWkf748Dj2Ty3JEmSGldEcN7onpw3uuc+r734k4eoXBHkfvbRJqjs4Pj8QUmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJElSIzhynq5nQJYkSVIjCBL1f4x1UzIgS5IkSRkMyJIkSVIGA7IkSZKUwYAsSZIkZTAgS5IkSRkMyJIkSVIGA7IkSZKyLlwHWZIkSfog10GWJEmSjkAGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJElZ9/46yCk1//WQDciSJElqFM0/GlczIEuSJEkZDMiSJElSBgOyJEmSlKHZB+SIGBgRP42I6RltfSPi0YiYFhE3NWV9kiRJOrpkNSDXBNi1EfHGXu2TI2JRRCw9UMBNKS1LKV21V/Ng4LGU0pXAsAYuW5IkSS1YtkeQ7wEmZzZERC5wB/BxqsPtpRExLCJGRsTv9/roWke/rwBTI+LPwP9msX5JkiS1MHnZ7DylNDMi+u/VPB5YmlJaBhARDwLnp5RuAc6tZ9dXAN+s6X868LO9d4iIa4FrAfr27Xtob0CSJEkN4v/WQYaIJi7mAJpiDnIvYEXGdklNW60ioigi7gLGRMQ/1DT/EfhiTfvy2o5LKd2dUhqbUhpbXFzcMJVLkiTpEAWJZp6Ma2R1BLkOtX1n6lw3OqW0Abhur7Y3gCkNXJckSZLUJCPIJUCfjO3ewKomqEOSJEnaR1ME5NnAoIgYEBH5wFTg0SaoQ5IkSdpHtpd5ewCYBQyJiJKIuCqlVAHcADwBLAAeTinNz2YdkiRJUn1lexWLS+tofxx4PJvnliRJkg5Fs3+SniRJktSYDMiSJEnKuj3rIDdxHfVhQJYkSVIjOTLWQTYgS5IkSRkMyJIkSVIGA7IkSZKUwYAsSZIkZTAgS5IkSRkMyJIkSVIGA7IkSZIaQc06yKn5r4RsQJYkSVLWBUfGQ0LAgCxJkiR9gAFZkiRJymBAliRJkjIYkCVJkqQMBmRJkiQpgwFZkiRJymBAliRJUiNIGf9t3gzIkiRJyrrqdZCjqcuoFwOyJEmSlMGALEmSJGUwIEuSJEkZDMiSJElSBgOyJEmSlMGALEmSJGUwIEuSJKnRpCNgIeRmH5Aj4oKI+HFE/DYizqppGxgRP42I6U1dnyRJkg7MdZBrRMS0iFgbEW/s1T45IhZFxNKIuGl/faSUfpNSugb4HHBJTduylNJVWStckiRJLVZelvu/B7gduO/9hojIBe4AzgRKgNkR8SiQC9yy1/FXppTW1nz9zzXHSZIkSVmT1YCcUpoZEf33ah4PLE0pLQOIiAeB81NKtwDn7t1HRATwXeAPKaWXs1mvJEmS1BRzkHsBKzK2S2ra6nIj8DFgSkRcBxARRRFxFzAmIv6htoMi4tqImBMRc9atW9dApUuSJOlol+0pFrWpbXZ2nfczppRuBW7dq20DcN3+TpJSuhu4G2Ds2LFHwP2SkiRJag6aYgS5BOiTsd0bWNUEdUiSJEn7aIqAPBsYFBEDIiIfmAo82gR1SJIkqdGkmv82/z/sZ3uZtweAWcCQiCiJiKtSShXADcATwALg4ZTS/GzWIUmSpKbX/KNxtWyvYnFpHe2PA49n89ySJEnSoWj2T9KTJEmSGpMBWZIkScpgQJYkSZIyGJAlSZKkDAZkSZIkKYMBWZIkSVkX76+DfASs9WZAliRJUqNIRFOXUC8GZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSlHWB6yBLkiRJRyQDsiRJkpTBgCxJkiRlMCBLkiRJGQzIkiRJUgYDsiRJkpTBgCxJkqRGkKr/m5q4jHowIEuSJCnrqtdBPjIYkCVJkqQMBmRJkiQpgwFZkiRJymBAliRJkjI0+4AcERdExI8j4rcRcVZGe9uImBsR5zZlfZIkSTq6ZDUgR8S0iFgbEW/s1T45IhZFxNKIuGl/faSUfpNSugb4HHBJxktfAx5u8KIlSZLUouVluf97gNuB+95viIhc4A7gTKAEmB0RjwK5wC17HX9lSmltzdf/XHMcEfEx4E2gIJvFS5IkqUYDLWCcjoDF3iJlebXmiOgP/D6lNKJm+0PAt1JKZ9ds/wNASmnvcPz+8QF8F3gypfRUTdvNQFtgGLAT+FRKqWqv464Frq3ZHAIsatA31nJ0AdY3dRHyOjQTXofmw2vRPHgdmgevw6Hrl1Iq3rsx2yPItekFrMjYLgEm7Gf/G4GPAR0j4riU0l0ppX8CiIjPAev3DscAKaW7gbsbrOoWKiLmpJTGNnUdLZ3XoXnwOjQfXovmwevQPHgdGl5TBOSopa3OYeyU0q3ArXW8dk8D1SRJkiQBTbOKRQnQJ2O7N7CqCeqQJEmS9tEUAXk2MCgiBkREPjAVeLQJ6lD9OE2lefA6NA9eh+bDa9E8eB2aB69DA8vqTXoR8QAwierJ42uAb6aUfhoR5wA/oHrlimkppZuzVoQkSZJ0ELK+ioUkSZJ0JGn2T9KTJEmSGpMBWZIkScpgQJYkSZIyGJB1UCJiYET8NCKm769NjScihkXEwxHxw4iY0tT1tFQRcWpE3BURP4mIF5q6npYqIiZFxLM112JSU9fTkkXE8TXXYXpEfL6p62mp/B19aAzILUhETIuItRHxxl7tkyNiUUQsjYib9tdHSmlZSumqA7WpfhrimgAfB25LKX0euDxrxR7FGuhn49mU0nXA74F7s1nv0aqBfh4SsA0ooHrdfR2CBvqZWFDzM/FXgE95OwTZ+r2tA3MVixYkIk6j+hfHfSmlETVtucBi4Eyqf5nMBi6legm+W/bq4sqU0tqa46anlD4wWllbm/avIa5JzedvAjuAk1NKpzRC6UeVBv7ZeBi4OqW0tZHKP2o00M/D+pRSVUR0A/4rpXRZY9V/NGmon4mIOA+4Cbg9pfSLxqr/aJHt39uqW1M8alpNJKU0MyL679U8HliaUloGEBEPAuenlG4Bzm3cClueBrwmX6j5R/ORbNV6NGuo6xARfYEthuND08D/Rm0CWmejzpagoa5FSulR4NGIeAwwIB8kf283HadYqBewImO7pKatVhFRFBF3AWMi4h/qatNhOdhr0j8i7gbuA/4zy7W1JAd1HWpcBfwsaxW1TAf783BhRPwIuB+4Pcu1tTQHey0mRcStNdfj8WwX14Ic9u9tHZgjyIpa2uqcd5NS2gBcd6A2HZaDvSbLgWuzVk3LdVDXASCl9M0s1dKSHezPwyP4l5RsOdhrMQOYka1iWrDD/r2tA3MEWSVAn4zt3sCqJqpF1bwmzYPXoXnwOjQfXovmwevQCAzImg0MiogBEZEPTAUebeKaWjqvSfPgdWgevA7Nh9eiefA6NAIDcgsSEQ8As4AhEVESEVellCqAG4AngAXAwyml+U1ZZ0viNWkevA7Ng9eh+fBaNA9eh6bjMm+SJElSBkeQJUmSpAwGZEmSJCmDAVmSJEnKYECWJEmSMhiQJUmSpAwGZEmSJCmDAVmSmkhE9I+IN5rgvGMj4tbGPm/NuR+PiGNqvt7WFDVI0oHkNXUBkqTGlVKaA8xponOf0xTnlaSD4QiyJB2kiPhNRMyNiPkRcW1G+7aIuDkiXo2IFyOiW037sTXbsyPi27WNnEZEbkT8Z80+r0XE3zRAnfdExF0R8WxELI6Ic2vaJ0XE7w+yr+KI+FVNfbMj4pSa9m9FxP0R8eeIWBIR19S094iImRExLyLeiIhTa9qXR0SXvfqOmvf+RkS8HhGXZNQ5IyKmR8TCiPh5RMThfl8k6UAcQZakg3dlSmljRLQBZkfEr1JKG4C2wIsppX+KiP/4/9u7n1CpyjiM49/HksQoF+bCRRZa4R+wq6Uo+W+RbVqYVqjLkkBBJSJaB9aqQC60ibIiQUGCEpK6Eqk3hbhRXl0IEtQmLhVEGeE1xZ4W5710kjNzPWNTBM8HBt4z5z2/933PYvjNO7+ZAZ4BXgIGgUHbByVt7xBzG3DB9jJJtwCnJB21/W29k6TPgNsarn/e9icNz98NrAXmAcck3dPDeilr2Gv7pKQ5VH9zu6CcWwysoFr/aUlHgK3AkO2XJd0ETO8SexMwANwP3EF1T4fLuSXAImAMOAU8BJzscQ0REdclCXJERHu7JW0s7TuBe4GfgMvAxM7sl8D60l4JPFbaB4BXG2I+AiyW9EQ5nlHi/i1Btr265VwP2f4D+FrSN8D8ltdPeBhYWNvAvV3SRKJ+2PY4MC7pGLAc+AJ4S9JU4APbo11irwIO2r4K/CDpBLAM+BUYsf0dgKRRqoQ/CXJE9FUS5IiIFiSto0oWV9q+KOk4MK2cvmLbpX2Vdq+xAnbZHppk/LY7yJ7kuB77baod27GGWuEpVGsev+aaxjFsD0taAzwK7Jf0iu13Ow3daU7A77V223saEdGT1CBHRLQzA/i5JMfzqUoLJvM58Hhpb+nQZwjYUXZckXSfpFuv7WR7te2BhkdTcgzwpKQpkuYBc4HznSZp+6kSq+mLdEeBnRMHkgZq5zZImiZpJrCOqkTiLuBH228A+4ClncYFhoHNpQ57FrAGGOnSPyKir5IgR0S08zFws6SzwB6q5HcyzwLPSRoBZgMXGvq8CZwDvio//fY6/8xu6XngBPARsN32pR7j7AYeLF8gPAfUa6lHgCNU92KP7TGqRHlU0mmqNweDXWK/D5wFzgCfAi/Y/r7HeUZE3DD99WlgRET0g6TpwLhtS9oCbLW94V8Y9x3gQ9vv9XGMF4HfbDfVVUdE/C+llisiov8eAF4rP1H2C/D0fzyfiIjoIjvIERERERE1qUGOiIiIiKhJghwRERERUZMEOSIiIiKiJglyRERERERNEuSIiIiIiJo/ARIN73MxBvIuAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 720x360 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plt.figure(figsize=(10, 5))\n",
"epsilons = np.logspace(-12, 0, 101)\n",
"ax = subplot(111)\n",
"errors = []\n",
"errors2 = []\n",
"axis = np.random.randn(3)\n",
"for epsilon in epsilons:\n",
" a = np.array([1., 1., 1., np.pi - epsilon])\n",
" #a = np.array([axis[0], axis[1], axis[2], np.pi - epsilon])\n",
" \n",
" a[:3] = a[:3] / np.linalg.norm(a[:3])\n",
"\n",
" R = matrix_from_axis_angle(a)\n",
"\n",
" angle1 = np.arccos((np.trace(R) - 1.0) / 2.0)\n",
"\n",
" r = np.array([R[2, 1] - R[1, 2], R[0, 2] - R[2, 0], R[1, 0] - R[0, 1]])\n",
" angle2 = np.pi - np.arcsin(np.linalg.norm(r) / 2.0)\n",
"\n",
" error = abs(a[3] - angle1)\n",
" error2 = abs(a[3] - angle2)\n",
" errors.append(error)\n",
" errors2.append(error2)\n",
"plot(epsilons, errors, label=\"via arccos\")\n",
"plot(epsilons, errors2, label=\"via arcsin\")\n",
"\n",
"ax.set_xlabel(\"angle = pi - epsilon\")\n",
"ax.set_ylabel(\"angle error [rad]\")\n",
"ax.set_xscale(\"log\")\n",
"ax.set_yscale(\"log\")\n",
"ax.set_ylim((10e-25, 10e0))\n",
"ax.legend()\n",
"plt.tight_layout()\n",
"plt.savefig(\"pi_angles.png\")\n",
"np.mean(errors2), np.median(errors2), np.max(errors2)"
]
},
{
"cell_type": "code",
"execution_count": 30,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(1.2834549483914063e-18, 0.0, 5.551115123125783e-17)"
]
},
"execution_count": 30,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAsgAAAFgCAYAAACmDI9oAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdd3yW1d3H8c/J3iEJATJJwgp7BRIZglgVtwVaAYvWRRER7aN9tH20jtbRVmvdiopbQVERBSeKOAhL9oYkkBAgIYHsfZ/nj4Q0QBICZBG+79crL7iv+1zX9bvvvgpfjr/rHGOtRUREREREKjm1dAEiIiIiIq2JArKIiIiISA0KyCIiIiIiNSggi4iIiIjUoIAsIiIiIlKDS0sX0Bzat29vo6KiWroMEREREWlFVq9efdBaG3zs8bMiIEdFRbFq1aqWLkNEREREWhFjzO7ajqvFQkRERESkBgVkEREREZEaFJBFRERERGo4K3qQRURERM4EZWVlpKWlUVxc3NKltCkeHh6Eh4fj6uraoPEKyCIiIiKtRFpaGr6+vkRFRWGMaely2gRrLVlZWaSlpREdHd2gc9RiISIiItJKFBcXExQUpHDciIwxBAUFndSsvAKyiIiISCuicNz4TvY7PSNbLIwx3sDzQCmwxFr7TguXJCIiIiJtRKuZQTbGzDbGZBhjNh5zfKwxZpsxZqcx5p6qw+OAedbam4Ermr1YERERkbNIeno6EyZMaOkymk2rCcjA68DYmgeMMc7Ac8DFQC9gkjGmFxAOpFYNq2jGGkVERETOOqGhocybN6/Rr1tR0TpjXKsJyNbapUD2MYeHAjuttUnW2lJgDnAlkEZlSIZW9BlEREREzmR33303zz//fPXrBx54gCeeeIKUlBT69OkDQEpKCiNHjmTQoEEMGjSIn3/+udZrXXXVVQwePJjevXsza9as6uM+Pj789a9/JT4+nmXLlrFy5UqGDRtG//79GTp0KHl5eVRUVHDXXXfRt29f+vXrxzPPPAPA4sWLGThwIH379uWGG26gpKQEgHvuuYdevXrRr18/7rrrrtP+Hlp7D3IY/50phspgHA88DTxrjLkU+LS2E40xU4GpAJGRkU1cpoiIiEjjevDTTWxOz23Ua/YK9eP+y3vX+f7EiRO54447mD59OgDvv/8+X3zxBQ6Ho3pMhw4d+Prrr/Hw8GDHjh1MmjSJVatWHXet2bNnExgYSFFREUOGDGH8+PEEBQVRUFBAnz59eOihhygtLSU2Npa5c+cyZMgQcnNz8fT0ZNasWSQnJ7NmzRpcXFzIzs6muLiY3//+9yxevJju3btz7bXX8sILL3Dttdfy8ccfs3XrVowxHD58+LS/p9Y++1rbI4fWWltgrb3eWntLXQ/oWWtnWWvjrLVxwcHBTVymiIiIyJlv4MCBZGRkkJ6ezrp16wgICDhuorGsrIybb76Zvn378pvf/IbNmzfXeq2nn36a/v37k5CQQGpqKjt27ADA2dmZ8ePHA7Bt2zZCQkIYMmQIAH5+fri4uPDNN98wbdo0XFwq53IDAwPZtm0b0dHRdO/eHYDrrruOpUuX4ufnh4eHBzfddBMfffQRXl5ep/09tPYZ5DQgosbrcCC9hWoRERERaTb1zfQ2pQkTJjBv3jz279/PxIkTj3v/ySefpGPHjqxbtw6Hw4GHh8dxY5YsWcI333zDsmXL8PLyYvTo0dXrEHt4eODs7AxUbuJR2xJstR231tZar4uLCytWrGDx4sXMmTOHZ599lm+//fakP3dNrX0GeSXQzRgTbYxxAyYCC1q4JhEREZE2a+LEicyZM4d58+bVunJFTk4OISEhODk58dZbb9X6oF1OTg4BAQF4eXmxdetWEhMTa71XbGws6enprFy5EoC8vDzKy8u58MILefHFFykvLwcgOzub2NhYUlJS2LlzJwBvvfUWo0aNIj8/n5ycHC655BL+85//sHbt2tP+DlpNQDbGvAcsA3oYY9KMMTdaa8uBGcCXwBbgfWvtppasU0RERKQt6927N3l5eYSFhRESEnLc+9OnT+eNN94gISGB7du34+3tfdyYsWPHUl5eTr9+/bjvvvtISEio9V5ubm7MnTuX2267jf79+3PBBRdQXFzMTTfdRGRkJP369aN///68++67eHh48Nprr/Gb3/yGvn374uTkxLRp08jLy+Oyyy6jX79+jBo1iieffPK0vwNT13R1WxIXF2drax4XERERaU22bNlCz549W7qMNqm279YYs9paG3fs2FYzgywiIiIi0hooIIuIiIiI1KCALCIiIiJSgwKyiIiIiEgNCsgiIiIiIjUoIIuIiIiI1KCALCIiIiL1Sk9Pr3XTkOY2bNiwZrmPArKIiIiI1Cs0NJR58+Y1+nVr24WvPj///HOj11AbBWQRERERAeDuu+/m+eefr379wAMP8MQTT5CSkkKfPn0ASElJYeTIkQwaNIhBgwbVGVqvuuoqBg8eTO/evZk1a1b1cR8fH/76178SHx/PsmXLWLlyJcOGDaN///4MHTqUvLw8Nm3axNChQxkwYAD9+vVjx44d1ecCLFmyhNGjRzNhwgRiY2O55ppraMzN71wa7UoiIiIi0ng+vwf2b2jca3bqCxc/VufbEydO5I477mD69OkAvP/++3zxxRc4HI7qMR06dODrr7/Gw8ODHTt2MGnSJGrbsXj27NkEBgZSVFTEkCFDGD9+PEFBQRQUFNCnTx8eeughSktLiY2NZe7cuQwZMoTc3Fw8PT158cUXuf3227nmmmsoLS2tdaZ5zZo1bNq0idDQUIYPH85PP/3EiBEjGuFLUkAWERERkSoDBw4kIyOD9PR0MjMzCQgIIDIykpSUlOoxZWVlzJgxg7Vr1+Ls7Mz27dtrvdbTTz/Nxx9/DEBqaio7duwgKCgIZ2dnxo8fD8C2bdsICQlhyJAhAPj5+QFwzjnn8PDDD5OWlsa4cePo1q3bcdcfOnQo4eHhAAwYMICUlBQFZBEREZE2rZ6Z3qY0YcIE5s2bx/79+5k4ceJx7z/55JN07NiRdevW4XA48PDwOG7MkiVL+Oabb1i2bBleXl6MHj2a4uJiADw8PHB2dgbAWosx5rjzJ0+eTHx8PAsXLuSiiy7ilVdeYcyYMUeNcXd3r/69s7Mz5eXlp/W5a1IPsoiIiIhUmzhxInPmzGHevHm1rlyRk5NDSEgITk5OvPXWW7W2P+Tk5BAQEICXlxdbt24lMTGx1nvFxsaSnp7OypUrAcjLy6O8vJykpCRiYmKYOXMmV1xxBevXr2/cD3kCCsgiIiIiUq13797k5eURFhZGSEjIce9Pnz6dN954g4SEBLZv3463t/dxY8aOHUt5eTn9+vXjvvvuIyEhodZ7ubm5MXfuXG677Tb69+/PBRdcQHFxMXPnzqVPnz4MGDCArVu3cu211zb656yPacwn/lqruLg4W1vzuIiIiEhrsmXLFnr27NnSZbRJtX23xpjV1tq4Y8dqBllEREREpAYFZBERERGRGhSQRURERFqRs6H9tbmd7HeqgCwiIiLSSnh4eJCVlaWQ3IistWRlZdW6HF1dtA6yiIiISCsRHh5OWloamZmZLV1Km+Lh4VG9qUhDKCCLiIiItBKurq5ER0e3dBlnPbVYiIiIiIjUoIAsIiIiIlKDArKIiIiISA0KyCIiIiIiNSggi4iIiIjUoIAsIiIiIlLDGRmQjTE9jTEvGmPmGWNuael6RERERKTtaPaAbIyZbYzJMMZsPOb4WGPMNmPMTmPMPfVdw1q7xVo7DfgtENeU9YqIiIjI2aUlZpBfB8bWPGCMcQaeAy4GegGTjDG9jDF9jTGfHfPToeqcK4AfgcXNW76IiIiItGXNvpOetXapMSbqmMNDgZ3W2iQAY8wc4Epr7aPAZXVcZwGwwBizEHj32PeNMVOBqQCRkZGNVr+IiIiItG2tZavpMCC1xus0IL6uwcaY0cA4wB1YVNsYa+0sYBZAXFycbaxCRURERKRtay0B2dRyrM5Qa61dAixpqmJERERE5OzVWlaxSAMiarwOB9JbqBYREREROYu1loC8EuhmjIk2xrgBE4EFLVyTiIiIiJyFWmKZt/eAZUAPY0yaMeZGa205MAP4EtgCvG+t3dTctYmIiIiItMQqFpPqOL6IOh64ExERERFpLq2lxUJEREREpFVQQBYRERERqUEBWURERESkBgVkEREREZEaFJBFRERERGpQQBYRERERqUEBWURERESkBgVkEREREZEaFJBFRERERGpQQBYRERERqUEBWURERESkBgVkEREREZEaFJBFRERERGpQQBYRERERqUEBWURERESkBgVkEREREZEaFJBFRERERGpQQBYRERERqUEBWURERESkBgVkEREREZEaFJBFRERERGpQQBYRERERqUEBWURERESkBgVkEREREZEaFJBFRERERGpQQBYRERERqcGlpQs4FcYYJ+BvgB+wylr7RguXJCIiIiJtRLPPIBtjZhtjMowxG485PtYYs80Ys9MYc88JLnMlEAaUAWlNVauIiIiInH1aYgb5deBZ4M0jB4wxzsBzwAVUBt6VxpgFgDPw6DHn3wD0AJZZa18yxswDFjdD3SIiIiJyFmj2gGytXWqMiTrm8FBgp7U2CcAYMwe40lr7KHDZsdcwxqQBpVUvK2q7jzFmKjAVIDIyslFqFxEREZG2r7U8pBcGpNZ4nVZ1rC4fARcZY54BltY2wFo7y1obZ62NCw4ObrxKRURERKRNay0P6Zlajtm6BltrC4Ebm64cERERETlbtZYZ5DQgosbrcCC9hWoRERERkbNYawnIK4FuxphoY4wbMBFY0MI1iYiIiMhZqCWWeXsPWAb0MMakGWNutNaWAzOAL4EtwPvW2k3NXZuIiIiISEusYjGpjuOLgEXNXI6IiIiIyFFaS4uFiIiIiEiroIAsIiIiIlKDArKIiIiISA0KyCIiIiIiNSggi4iIiIjUoIAsIiIiIlKDArKIiIiISA0KyCIiIiIiNSggi4iIiIjUoIAsIiIiIlKDArKIiIiISA0KyCIiIiIiNSggi4iIiIjUoIAsIiIiIlKDS0sXICICsPdwEcuTsght50lCTFBLlyMiImcxBWQRaXbWWtIOFZGYlMXy5GwSk7JIO1QEQJC3G6vu/RXGmBauUkREzlYKyCLS5Ky17M4qZHlyFsuTKgNxek4xAAFergyNDuSG4dEcLizl6W93sjMjn24dfVu4ahEROVspIItIo7PWknywgMSk7OpQvD+3MhAHebsRHxPIH6KDSIgJolsHH5ycKmeL92QV8vS3O0lMylJAFhGRFqOALCKnzVrLrsx8Eqtmh5cnZ5OZVwJAex93EmICiY8JIiE6kK4dfOpsn4gI9CTU34PEpGymnBPVjJ9ARETkvxSQReSkORyWHRn5LE/OIjEpixXJ2RzMLwWgo587w7oEER8dRHxMIDHtvRvcT2yMISEmiO+3Z2KtVR+yiIi0CAVkETkhh8OydX9edbvE8uQsDhWWARDq78HIbsHERweSEBNE5yCv0wq2CTFBfLRmr/qQRUSkxdQZkI0x4xpwfrG1dlEj1iMirUCFw7JlX251u8SK5GxyiioDcXiAJ2NiOxIfE8g5MUGEB3g26kzvkSXe1IcsIiItpb4Z5JeBT4D6/uY7F1BAFjnDlVc42JSeWz1DvCIlm7zicgA6B3lxUe+O1S0T4QFeTVpLdR9ysvqQRUSkZdQXkD+31t5Q38nGmLcbuR4RaQZlFQ427s2pXmViVcoh8ksqA3FMe28u6xdSHYhD/D2btbYjfchLd6gPWUREWkadAdla+7sTndyQMSLS8krLHWzYe7h6lYnVuw9RWFoBQJdgb64cEFq9ykQHP48WrhbiYwL5aM1edmXm07VD3W0W2QWleLo64+nm3IzViYhIW3fKPcjW2o8avxwRaQwl5RWsT8shcVdlD/Hq3YcoKqsMxN07+jBhcDjx0UEMjQ4k2Ne9has93pE+5GVJ2XUG5JLyCi556gcGdw7guWsGNWd5IiLSxtXXYnF51a8dgGHAt1WvzwOWAArIIq1EcVkFa1MPV68wsXr3IUrKHQDEdvLl6iERJMQEMiQqkCCf1heIjxUZ6EWIvweJSVlMSehc65gFa9PZn1vMoo37SDlYQFR772auUkRE2qr6WiyuBzDGfAb0stbuq3odAjzXPOXVzhjTC3gAyAIWW2vntWQ9Is2tuKyCX3YfIjE5m+VJWaxJPUxpuQNjoFeIH9fEdyYhJpCh0YG083Jr6XJP2pE+5B/q6EO21vLqj8lEBXmRfriYV39M5m9X9WmhakVEpK1pyDrIUUfCcZUDQPdTvaExZjZwGZBhre1T4/hY4CnAGXjFWvtYPZe5GHjGWvuDMWYBoIAsbVphaTm/7D5ctexaFutScyitcOBkoE+YP9ed05n46CCGRAXi7+Xa0uU2ioSYQD6uow952a4stu7P45/j+7F69yE+WJ3KHy/oTqB34/1joKCknE/WpvObuHBcnZ0a7boiItL6NSQgLzHGfAm8B1hgIvDdadzzdeBZ4M0jB4wxzlTOSl8ApAErq4KvM/DoMeffALwF3G+MuQIIOo1aRFqlgpJyVu0+VBmIk7JYn5ZDucPi7GToE+bP9SOiSIgOYnBUAH4ebSMQH6u+PuRXfkymvY8bVwwIZWBkO+auSuXtxN3MPL9bo93/ue928vySXRSXVXDDiOhGu66IiLR+JwzI1toZVQ/sjaw6NMta+/Gp3tBau9QYE3XM4aHATmttEoAxZg5wpbX2USpnm2tza1WwrrUX2hgzFZgKEBkZearlijSLvOIyVqUcIjE5i8SkbDbuzaHCYXFxMvQL9+fmc2NIiAlicOcAfNzPjg0w6+pD3pWZz7dbM7j9/G54uDrTraMvY2I78MbPKUw9NwYP17pXtNi4N4fuHX1xc6l/Rjgrv4TXf07BGHj62x2MHxyOv2fb/IeIiIgcr0F/01atWNGUD+WFAak1XqcB8XUNrgrYfwG8gX/VNsZaOwuYBRAXF2cbqU6RRpFTVMbK5MoH6pYnVwZihwVXZ8OAiHZMH92F+OggBnVuh5fb2RGIj1VXH/JrPyXj5uzE72qE5ptHxjDp5UQ++mUvk+Nr/wfx3JV7uPvDDUyOj+SRX/et994vLU2iuKyCpycOZOacNTz/3U7+fEnPxvtwIiLSqp3wb15jTALwDNATcKOy7aHAWuvXiHXUthNAnaHWWptC1eywyJngcGEpy5Ozq1eZ2LwvF2vBzcWJgRHtmDGmGwkxgQyKDKh3BvRsc2wf8uHCUuatTuOqgaFHLU+XEBNI3zB/XvkhiYlDInByOvqPlBXJ2dw7fyN+Hi68u3wP4weFMbhzYK33zMgr5s1lKVw1MIzL+4fy3bYMXvs5hSnndG7yXQRFRKR1aMjU1LNU9h1/AMQB1wJdG7mONCCixutwIL2R7yHSbLLyS1iRnM3y5MqNObYdyMNacHdxYlBkALef342EmCAGRLRTIK7HsX3I7yzfQ3GZ47ieYGMMU8+N4bb31rB4awYX9OpY/V5qdiHT3l5NRIAX79wcz4QXlvHnjzbw2W0ja221eGHJLsoqLDPHVPYz33VhDxau38fjX27jPxMHNuGnFRGR1qKhLRY7jTHO1toK4DVjzM+NXMdKoJsxJhrYS2Ugn9zI9xBpMpl5lYH4yCoT2w/kA+Dh6kRc50Au7RtCfEwQ/SP8cXdRIG6oI33Iy5OyuDougjeXpTCia3tiOx3/H7Au7tOJsHaevLw0qTog55eUc9MbqyivcPDKdXGE+Hvy0JW9ufGNVbz8QxK3nnf0v/X35xTzTtUM85F1lUPbeXLDiGheWLKLG0fE0Dfcv8k/t4iItKyGBORCY4wbsNYY809gH5W9v6fEGPMeMBpob4xJA+631r5qjJkBfEllC8dsa+2mU72HSFPLyC2uXoM4MSmLXZkFAHi5ORMXFciVA8Kq/rN/uxM+ECZ1+28f8kEWbdjHgdwSHhvXr9axLs5O3Dgimoc+28yaPYfoF96OO+asYWdmPq9fP4SYYB8Azu/ZkYv7dOLpxTu4rF8InYP++8fZc9/txOGw3Dbm6NUwbhndhbkrU3lk0RbevTn+uHWZ61NcVkFBSfkZsUGLiIhUakhAngI4ATOAP1LZCjH+VG9orZ1Ux/FFwKJTva5IU9qXU1TdP7w8KZukg5WB2MfdhbioACYMrtyprk+Yv9bMbWTx0ZV9yI99vpUuwd6M6h5c59jfDongP99s5+Ufkugc5M03WzJ48IrejOx29Dn3X96bH3Yc5N75G3nzhqEYY9h7uIg5K/fw2yERRAQe3Wvs5+HKzDFdeeDTzXy3LYMxsR1piKLSCn770jJ2ZxWwYMaIRt/tz1rL+6tSGdal/XE1i4jIqas3IFcto/awtfZ3QDHwYLNUJdLC9h4uInFXVvUqE7uzCgHw9XBhaFQgE4dGkBATRK8QP1wUiJvUkT7k/bnFPPzrPsc9gFeTj7sL1yR05sXvd2EtTI6P5Npzjt+qupO/B3+6qAf3L9jEgnXpXDkgjGe/3YnBMOO82h+xmBzfmdd/TuHRRVs5t1vwCf93t9byvx+uZ2N6Dj5uLkx9axUfTx+OdyMu07dgXTp3f7iBqCAvPp4+nIBG3ChFRORsVu+f1NbaCmNMsDHGzVpb2lxFiTQnay1ph4pITKpcg3h5chZph4oA8Pd0ZWh0IFMSOpMQE0TPED+c6wlo0vg6B3nRyc+D4vIKxg0MP+H43w+LYvaPyQyMbMeDV/Susx3idwmd+WjNXh76dDMx7X34YFUq18RHEtrOs9bxbi5O3D02llve+YUPVqcxaWj966u/8P0uPl2Xzp8u6kH/8HZcO3s5f5q3jucmDzqpFo26FJaW8+iirdXbbd/yzmrevCFeLT0iIo2gIVMZKcBPVTvbFRw5aK39d1MVJdKUrLXsziqsbpdITMoiPacYgEBvN4ZGBXLjiGjio4OI7eRb74ylND1jDA9c0RsXJ4On24kfcOzo58HiO0cR7Oteb7uLs5PhkV/34Ypnf2LirGU4ORmm1zF7fMTYPp0Y3DmAJ77aRs8QPwZEtKt13OItB/jXl9u4vH8o00d3wRjDPRfH8siirTy/ZNdxDweeiue/28X+3GLmTTuHtENF3DF3LffO38A/xvdrlAAuInI2a0hATq/6cQJ8TzBWpNWx1pJ0sKC6hzgxKYsDuSUAtPdxIz46iGkxgSTEBNE12EeBuBUa26fTSY1v6HrFvUP9uXFENLOWJnHD8Gg6+nnUO94Yw9+u7MMNr69k3PM/cf3waO68sPtRm7nszMjj9jlr6R3qxz9rhNWbR8awYW8uj3+1jV6hfpzXo8NJfaaa9mQVMuuHJK4cEEpcVCBxUZU7DD7z7U66dvBh6rldTvnaIiICxtq2v8lcXFycXbVqVUuXIc3EWsvOjPzqVSaWJ2eTmVcZiIN93UmICSI+OpCEmEC6BPtotu0sV1RawZyVexg/OBw/j4ZtJ51XXMY/v9jGW4m7CQ/w5OFf92VU92ByCsu48rkfyS8pZ8GMEce1axSWljPu+Z9JP1x0Wg/t/eGtVfyw4yDf3jmaTv6Vod7hsNz23hoWbdzHrClxR60FLSIitTPGrLbWxh13vK6AbIx5wFr7wAkuesIxrYECctvmcFi2Z+QdtcpEVkFly3wnPw/iYwKJjw4iPiaQmPbeCsTSaFamZHP3h+tJyizg1wPDyMwrYXlyFu/dnEBcVO079aVmF3L5sz/Swdf9lB7a+3HHQX736nL+dFGP41o1ikoruHrWMnZm5PPBtHPoHXrqazZba3n8q20sXL+PR8b1ZViX9qd8LRGR1upUAnIaUF+fsQFuttbGNk6JTUcBuW1xOCxb9+dVb8qxIjmbQ4VlAIS18yQ+OpD4qpaJyEAvBWJpUsVlFTz33U5eWLKLcoflH+P7cvWQ+h/g+3HHQa6dvZxOfh6EB3gR6O1GoI8bQd5uBHpXtv30Cj1+M5SyCgeXPPUDJeUOvvrjubXuwpiRW8yVz/0EwKPj+jKia/uTXmmlwmG5d/5G3luxB39PV/KKy5h5fjduG9NND6mKSJtyKgH5/gZcN99a+8TpFtfUFJDPbBUOy5Z9udWrTKxMySanqDIQRwR6Vs4OR1cGYq0FKy1l+4E8dmXkc3HfkAaN/2x9Oos27CMrv5TsgsqfQ4WlOKr+SL6wV0du/1W3o2aBX/spmQc/3cysKYO5sHfdfdmb0nOY8uoKsgtKCfJ249J+IVzRP5RBkQEn7LEvLXfwP++v5bP1+5g+ugu3nteV+z7ZyEe/7CUhJpCnJg48Ya+2iMiZ4qQDcluigHxmKa9wsCk9t7pdYkVKNnnF5UDlkl8JVe0S8TFBhNWxJJfImajCYcnKL+G9Fam88mMSecXl1UG5k58Hox9fwoCIdtWbm9SnpLyCJdsyWbAunW82H6Ck3EFYO08u7x/KFf1D6Rnie9w1ikoruOWd1SzZlsk9F8cybdR/H/abtzqN++ZvxMvNmSevHsC59WzYIiJyplBAVkButcoqHGzcm1O9BvGqlEPkl1QG4pj23tXtEvHRQdUPJIm0dTlFZbz+U0p1UA4P8GR/TjFf3DGSrh1ObkGh/JJyvt68n0/WpvPDjoNUOCxdgr25vH8ol/cPpUuwD7nFZdz0+ipW7s7m4av6Mjn++DaRHQfymPHuGrYdyOP3w6L4XULkSdciItKaKCArILcapeUONuw9TGLVGsSrdx+isLQCgC7B3pVhOCaIhOhAOug/5cpZLqeojNd+SubVH5OZktCZ/x17eo99ZBeU8vnGfSxYm86KlGyshd6hfpRXWHZl5vPk1QO4vH9onecXlVbwt4WbmbNiDw4LfcL8uGpAGFf0D9X/X0XkjHNKAblqq+mZ1tonm7K4pqaA3LJKyitYl5rD8qQsEpMrA3FxmQOA7h19qmeHh0YHEuzr3sLVirROR/6sbsyHTvfnFLNwwz4+XZfOjgN5PDt5EOfFNmx95oy8Yj5dt49P1u5lfVoOTgaGdWnP+T070K2DL106eNPJz0MPyYpIq3bKM8jGmCXW2tFNVVhzUEBuXsVlFaxNPVy5ykRSNo7aW+4AACAASURBVL/sOURJeWUgju3kW70O8dDoQIJ8FIhFWgOHw57yJjk7M/L5ZO1e5q/dS2p2UfVxbzdnunTwoUuwDyO7teeK/qEnvaKGiEhTOp2A/DDgD8zl6K2mf2nsIpuKAnLTKiqtYM2eQ5WrTCRnszb1MKXlDoyBXiF+RwXidl5uLV2uiDQRay2ZeSXszMxnV2YBuzLy2ZWZz/YDeRzILSG6vTczz+/KFf3DtFyciLQKpxOQv6vlsLXWjmms4pqaAnLjKiwtZ/XuQyyv6iFel3aYsgqLk4E+Yf7VgTguKhB/z4btTCYibZe1li83HeA/32xn6/48ugR7M/P8blzWL1RBWURalB7SU0A+Zfkl5axKyWZ51dbN69NyKHdYnJ0MfcP8q1eZiOscgG8Dt+oVkbOPw2H5YtN+/vPNdrYfyKdbBx8eG9+XwZ1r33XwWGUVDlzVoiEijeh0ZpD9gfuBc6sOfQ88ZK3NafQqm4gC8snJLS6rDMRJ2SQmZ7Nxbw4VDours6FfeLuqneqCGNw5AJ+T3CZXRMThsCzauI/HPt9KYWkFC2eOIMS//jXNVyRnc/1rK7jvsl5MHFr/ToUiIg11OgH5Q2Aj8EbVoSlAf2vtuEavsokoINcvp7CMFSmVs8PLk7PZlJ6Dw4Krs2FARLvqVSYGdW6Hl5sCsYg0jl2Z+VzxzI/0CvXjvZsT6nyAb19OEZc/8yMH80vxcXfhm/8Z1aA10fOKy/RftUSkXqcTkNdaawec6FhrpoB8tEMFpVWBuLKHeMv+XKwFNxcnBka0q16DeGBkAJ5uzi1droi0YZ+s3cvtc9YybVQX7rn4+DWei8squPqlZezKLODpSQO45e1fOK9HB16cMrje6774/S7+8cVW/nxxLDePjDnhcnMrkrPZfiCPa+IjtTSdyFmkroDckOnAImPMCGvtj1UXGg4UneAcaUWy8ktYkVzZQ5yYlMXW/XkAuLs4MSgygDvO7058TCADItrh4apALCLN58oBYSQmZfPi97sYGh3AmNiO1e9Za/m/jzeyLi2HWVMGMya2I3f8qjv/+GIrX23az4W9O9V6zSXbMvjHF1tp7+POI4u2kppdxP2X96p1htrhsLy4dBePf7kNh63c6v73w6Ob7POKyJmhIQF5GvBmVS8ywCHguqYrSU5XZl4Jy5OzqmeId2TkA+Dp6szgzgHceUEICV2C6Bfuj7uLArGItKz7L+/F2tTD/M/761g4cyRh7Sr7kd/4OYUPf0nj9vO7VYfhm0ZG88navdy/YBPDurY/7jmI3VkFzHxvDT06+vLhLcN4evEOXlqaRPrhIp6eNBDvGuNzCsu484O1fLMlg8v6hVTtEriF7p18GdalffN9ASLS6pxoJz0nYIK19n1jjB+AtTa3uYprLG29xeJAbnHlphxVq0zsyqxcrtrbzZnBUYHERweSEBNI37B2uLnoCXARaX2SDxZw+TM/0r2jD3P/cA6rUg7xu1eXc16PDsyaMvioTUzW7DnEuBd+5rpzonjgit7VxwtLyxn3/M/syynm0xkjiAzyAuCtxN3c/8lGeoX6Mfu6IXTw82Dj3hxueWc1+w4Xc++lPbluWBT5JeVc9dxPZBeUsmDGCCICvZr9exCR5nU6PchLrbXn1juolWtrAXlfTlH17PDy5GySD1YGYh93F+KiAqrXIe4T5q8lkUTkjPHZ+nRmvLuGCYPD+W5rBu28XJl/6/BaH7S7/5ONvJm4m4+nD2dARDustcycs5bP1qfz+vVDGdU9+Kjx3249wIx31xDg5cbk+EieWryDQC83nrtmEIM7B1SPS8rM58rnfiI8wIuPbhmm5zBE2rjTCcj3UdlzfOxOetmNXWRTOdMDctqhwqMC8Z7sQgB8PVwql1yLDiI+JpBeIX7axlVEzmj3zd/IW4m78XV3Yf6M4XQJ9ql1XF5xGRf8eykB3m4smDGc139K4eFFW/jfsT2YPrprreds3JvD9a+vJDOvhOFdg3hq4kDa17Ld/XdbM7jhjZVc2jeEZyYN1EN7Im3Y6QTk5FoOW2ttTGMV19TOpIBsrSU1u4jEGj3Eew9XPhPp7+nK0OjA6hniniF+2oVKRNqUkvIKHvp0M5f0DWF41/r7gL/YuJ9pb6/msn4hLNqwj7F9OvHc5EH1Btp9OUUkJmWdcLvr55fs5J9fbOPusbHcMrrLKX8eEWndTikgV/Ugn2Ot/akpi2tqrTkgW2tJySqsXoM4MSmLfTnFAAR6uzE0qrJ/OD4miB4dfY/qwxMROdtNfXMVX20+QPeOPnw8ffhRD+GdDmstM95bw6IN+3jk133p2sEHDxdnPFyd8HB1xt3VCYOhtMJBWbmDsgoHJVW/dvL3OOHGJyLSOpzSMm/WWocx5nHgnCarrAGMMTHA/wH+1toJVceuAi4FOgDPWWu/asESG8xay67MgqNWmcjIKwGgvY8b8dFB1YG4a7CPArGISD3+dlUffD1cuW1M10YLxwDGGP41oR9JmQX8+aMNJ31+Jz8PBkS0Y2BkOwZEtKNvuD/OTobMvBIO5JaQmVdMRl4JuUVljBsUTmg7BWqR1qQhLRYPAuuBj+yJBtd+/mzgMiDDWtunxvGxwFOAM/CKtfaxBlxr3pGAXONYAPC4tfbGus5ryRlkay07M/JJTMoiMblyc46D+ZWBuIOvO/FV7RIJMYF0CfZRr5uISCtSVFrBhr05FJdVUFLuoLisovKn3AGAm7PBzcUJV2cn3Jwrf92dVcCa1MOsTT3M7qzKZ0acDDjq+Bs0ItCTD/4wrEG7A4pI4zqdHuQ8wBuooPJhPUNlD7JfA298LpAPvHkkIBtjnIHtwAVAGrASmERlWH70mEvcYK3NqDqvtoD8BPCOtfaXumpo7oBcXuHgneV7SEzKYkVyNlkFpUDljEJ8zH97iKPbeysQi4i0YVn5JaxLO8y61BxcnAwd/Nzp4OtR/WvaoUKmvLqCDn7uzJ16DsG+xz80KCJN55QDciPdPAr4rEZAPgd4wFp7UdXrPwNYa48Nx8deZ16NFgsDPAZ8ba39ppaxU4GpAJGRkYN3797daJ/nRKy1nPPotzg7mcpAXLXKRGSglwKxiIgcZWVKNte+uoLOQV7MmZpAOy+3li5J5KxRV0A+4ZpgptLvqpZ7wxgTYYwZepr1hAGpNV6nVR2rq4YgY8yLwMAjYRq4DfgVMMEYM+3Yc6y1s6y1cdbauODg4GPfblLGGL64YyQ/3TOGf/92AL8dEkHnIM0Wi4jI8YZEBfLytXEkHSzg2tkryC0uO27Mzow87p2/gTGPL2HJtowWqFLk7NKQRXOfp/IhvclVr/OB507zvrUlxTqnsq21WdbaadbaLkdmma21T1trB1cdf/E062l0mgEQEZGGGtGtPS/+bhBb9uVy/WsrKSgpx+GwLN5ygCmvLudX/17K+6vSKCl3cMPrK3nj55SWLlmkTWvII7/x1tpBxpg1ANbaQ8aY001/aUBEjdfhQPppXlNEROSMNSa2I09NHMiMd39h8suJHC4qY3dWIZ38PPjTRT2YOCQCD1dnbp+zlvsXbGJXZj5/vayXNogSaQINCchlVQ/VWQBjTDDgOM37rgS6GWOigb3ARP47Qy0iInJWuqRvCE/8tj93fbCegRHt+NNFPbiodydca4Tgl6YM5h9fbGXW0iR2ZxXyzOSB+NWyHbeInLqGrGJxDXA1MAh4A5gA3Gut/aBBNzDmPWA00B44ANxvrX3VGHMJ8B8qV66Yba19+FQ/xIm05o1CREREjlVYWo6XW/1zWHNW7OHe+RuJCfbm1euGEODtxsG8Eg7mV/5k5pXQ0c+DC3t3aqaqRc48p7WKhTEmFjifyt7hxdbaLY1fYtNRQBYRkbbo550Hmfb2anKLy+sc868J/fhNXESd74uczU5pJ70jrLVbga2NXpWIiIicsmFd2zP/1uF8+Esavh6utPdxJ9jXnfY+bgR5u3PnB2v5v/kb6dbRlwER7Vq6XJEzRrOsg9zSNIMsIiJno+yCUi5/5kcqHJYFtw2ng6926xOp6ZTXQRYREZEzU6C3G7OuHczholJufecXSstP7hl7ay0PLNjEu8v3NFGFIq1TgwKyMaazMeZXVb/3NMb4Nm1ZIiIi0hh6h/rzj/H9WJlyiL99tvmkzp2zMpXXf07h319vp7zidBewEjlzNGQnvZuBecBLVYfCgflNWZSIiIg0nisHhDH13BjeStzN3JUNmw1OPljAQ59uppOfBwfzS1i6I7OJqxRpPRoyg3wrMBzIBbDW7gA6NGVRIiIi0rj+96IejOjanvvmb+KXPYfqHVtW4eCOuWtxc3Hig2nnEODlyoe/7G2mSkVaXkMCcom1tvTIC2OMC/VsCy0iIiKtj4uzE89MGkhHf3dufH0lP+88WOfYZ77dybrUwzw6ri8RgV5cOSCMrzcfIKewrBkrFmk5DQnI3xtj/gJ4GmMuAD4APm3askRERKSxBXi78eYN8QT5uDNl9gpe+SGJY1ezWr37EM9+u4Nxg8K4pG8IAOMHhVNa7uCzDektUbZIs2tIQL4HyAQ2AH8AFgH3NmVRIiIi0jSi23sz/9bh/KpnB/6+cAt3zF1LUWkFAPkl5fxx7lpC23ny4BW9q8/pE+ZHtw4+fLg6rd5rVzgsL36/i5SDBadVY3FZBat3Z1NcVnFa1xE5VSfcKMRa6wBervoRERGRM5yPuwsvXDOYF77fxeNfbWPHgXxemjKYZ7/dSdqhQuZMPQdfD9fq8cYYxg8O57HPt5KUmU9MsE+t1313xR4e+3wr327JYO4fEjDGnHRtaYcKmfb2ajbuzcXbzZkLe3fi8v4hjOgajJuLVqeV5lHnRiHGmA3U02tsre3XVEU1Nm0UIiIiUrvvtmYwc84arK2cQb71vC786aLY48YdyC3mnEcXM310V+66qMdx7x/ML2HM40twcjIcLixj9u/jGBPb8aRq+XnnQWa8t4aycgd3XdSDLfty+XzjfnKKyvD3dOXiPp04L7YDQd5u+Hq44uvhgq+HC95ulfN92YWlZOaVkJlXQkbVrz1DfBndQ2sLSO3q2iikvoDcub4LWmt3N1JtTU4BWUREpG4pBwuY9vZqvNycmTP1nDpnaq+dvYJdGfn88L/n4eR09OzwXR+s45O1e1kwYwTT3l6Np6szC2eOxNnpxLPI1lpe/TGZRxZtISbYh1lTBlfPUpeWO/hxZyafrtvHV5v2U1B6fNuFMWAARy2RxsPVicV3jiasneeJvwg569QVkOtssTiTArCIiIicuqj23nx++0gclnoD7fhBYdw+Zy2JyVkM69K++viqlGzmrU5j2qgu9Azx484LezDzvTV8snYv4waF13vvotIK7v5wPQvWpTO2dyce/21/fNz/G0/cXJwYE9uRMbEdKS6rYMu+XHKLy8kvLievuIy8ql8tEOzrTrCPO8G+7nTw9aDM4eDSp3/g0UVbeHbyoNP+nuTsccIeZGNMHse3WuQAq4A7rbVJTVGYiIiINB9jDM4nmOy9qHcnfN1d+HD13uqAXF7h4N75Gwn192Dm+V0BuKxvCLOW7uKJr7Zzab8Q3F2ca73e/pxirn99JVv35/Kni3owfXSXevuWPVydGRgZcFKfa9qoLvznmx1ce042Q6MDT+pcOXs1pNv938CfgDAqd9G7i8oH9uYAs5uuNBEREWlNPFydubRfCJ9v3EdBSTkAbyzbzdb9efz18l54VfUCOzkZ7h4by97DRbydWPvOfRm5xUx6OZHU7EJe+/0Qbj2v6yk91Hcifzi3C6H+Hjz46SYqauvBEKlFQwLyWGvtS9baPGttrrV2FnCJtXYucHL/jBMREZEz2vjB4RSWVvDFxv0cyC3mya+3M6p7MBf17nTUuJHdghnRtT3PfruD3OKjNxjJyKsMxxm5xbxxw5AmfYjO082ZP1/Sk03puXywKrXJ7iNtS0MCssMY81tjjFPVz29rvKd/iomIiJxF4joHEBnoxYe/pPHwwi2UVjh48Iretc7+3j02lkOFZby89L/dmAfzS7jm5eWkHy7mteuHMrhz07c9XNYvhCFRAfzry23HhXWR2jQkIF8DTAEygANVv/+dMcYTmNGEtYmIiEgrY4xh3KAwft6VxYJ16Uwb1YWo9t61ju0b7s9l/UJ45YdkMvKKyS4o5ZqXl5N6qJDZvx/SbD3Bxhjuv7w32YWlPLN4R7PcU85sJwzI1toka+3l1tr21trgqt/vtNYWWWt/bI4iRUREpPUYX7UyRWSgF9NHd6l37F0X9qCswsEjC7dwzSvLSckq4NXrhnBOl6DmKLVanzB/ro6L4LWfUtiVmd+s95YzT0NWsQgGbgaiao631t7QdGWJiIhIaxUR6MXfr+pDv3B/PFxrX6HiiKj23kwaGslbibtxc3HilWvjGN61fb3nNJU7L+zBwvX7+Ptnm3nt+qEtUoOcGU4YkIFPgB+AbwBtii4iIiL8LqHe/cSOMvP8buzOLuSG4VGc2z24CauqX7CvOzPP78bDi7bw3dYMzottmocDDxWUEuDt1iTXluZR50561QOMWWutHdBM9TQJ7aQnIiIiULkz39inlpJTWMacqQl06+h7wnPSDhXSyc8DF+cTP7r16bp0Zs5Zw+MT+jN+cP2bpEjLq2snvYY8pPeZMeaSJqhJREREpFm5uTjx6nVDcHIyTH5lOUn19CNba3ntp2TO/ed33D53LSeaVMzMK+Gvn2zEWvj7ws1kF5Q2dvnSTBoSkG+nMiQXGWNyjTF5xpjcpi5MREREpClEt/fm3ZvicTgsk19ezu6sguPGlFU4uO+TjTz46Wai2nuzcP0+3lle+6YnUBmm75u/kYKSCp6bPIi84nIeWbSlKT+GNKGGrGLha611stZ6Wmv9ql77NUdxIiIiIk2hW0df3r4pnuLyCia/vJy0Q4XV7+UUlXHD6yt5O3EP00Z14as7zmVU92Ae+mwzm9Jzar3ewg37+GLTfu64oBuX9gvhD6NimLc6jWW7sk6pvsLSctamHiYjtxiHdgBsdifsQQYwxgQA3QCPI8estUubsK5GpR5kERERqc3GvTlMfjmRdl5uzP1DAqXlDm54fSV7sgt5+Nd9+W1cBABZ+SVc8vQPeLu5sOC2Efi4/3edg6z8Ei54cinhAZ58dMswXJydKC6r4MInl+LiZPj8jpG4u9S/2kdN1lpufGMV327NAMDdxYmwAE/CA7wID/Dkot6dGNWCDzu2Jafcg2yMuQlYCnwJPFj16wONXWA9948xxrxqjJlX49hoY8wPxpgXjTGjm6sWERERaVv6hPnz5o3xZBeUMnFWIlc99xNZBaW8dWN8dTgGCPJx56mJA0nJKuDejzcc1Y/81wWbyC8u518T+lc/yOfh6szfr+pD0sECXliy66Rq+mB1Gt9uzWDquTH87cre/H5YFLGdfDlcWMqn69K5+c1VHMgtbpwvQGrV0B7kIcBua+15wEAgsyEXN8bMNsZkGGM2HnN8rDFmmzFmpzHmnvquUbVRyY3HHgbyqZzRTmtILSIiIiK1GRDRjtevH0JmXgkB3m7Mnz6chJjjNzJJiAnijl91Z/7adD5YVRk/Fm3Yx8L1+7j9V93o0enoFTHO7R7MlQNCef67XQ3enCT9cBF/+3Qz8dGB3DM2linnRPHnS3ry/DWDWTBjBAtvG4nDYU86dMvJaUhALrbWFgMYY9yttVuBHg28/uvA2JoHjDHOwHPAxUAvYJIxppcxpq8x5rNjfupaoPAHa+3FwN1UzmqLiIiInLK4qECW/Gk0C28bWefW2QC3nteV4V2D+OuCjSQmZXHf/I30DfPnD+fG1Dr+3kt74eHqxP8dM+tcG2std3+4ngpr+deE/jg5mePGRAZ5MX5QOO+u2MP+HM0iN5WGBOQ0Y0w7YD7wtTHmEyC9IRev6lPOPubwUGBn1cxwKTAHuNJau8Fae9kxPxl1XNdR9dtDgHtDahERERGpTwdfDzzd6u8VdnYyPHn1AHzcXZj0ciK5xWX86zf96lwjOdjXnXsu7kliUjYf/rK33mu/u2IPP+w4yF8u6UlkkFed42aM6YrDYXnxe80iN5WGrGLxa2vtYWvtA8B9wKvAVadxzzAgtcbrtKpjtTLGBBljXgQGGmP+XHVsnDHmJeAt4Nk6zptqjFlljFmVmdmgjhARERGRE+rg68F/rh6IkzHc8avuxHaqf3GviUMiGNw5gIcXbmbp9tozSWp2IQ8v3MKIru25Jj6y3utFBHoxYbBmkZtSg1axOK0bGBMFfGat7VP1+jfARdbam6peTwGGWmtva6oatIqFiIiINLacwjL8vVwbNHZnRj7Xv76C1OwiRnZrz18u6UnPkMpg7XBYJr+SyMa9uXz5x3MJa+d5wuulZhdy3uNLuCY+kgev7HNan+Nsdjo76TW2NCCixutwGtiyISIiItJaNDQcA3Tt4MM3/zOKey/tyfq0HC55+gfu+mAd+3KKeCtxN4lJ2dx3Wc8GhWOonEX+TVw4761IZV9O0al+BKlDSwTklUA3Y0y0McYNmAgsaIE6RERERJqNu4szN42MYemfzuPmkTEsWJvOeY8v4ZFFWxjdI/ioZeUaYvrorjis5fnv1Ivc2Jo0IBtj3gOWAT2MMWnGmButteXADCrXU94CvG+t3dSUdYiIiIi0Fv5ervzlkp4svnMUF/bqRCd/Dx4b1w9jjl+1oj6Vs8gRzF2ZSvphzSI3pibvQW4N1IMsIiIibVHaocpe5KuHRPD3q/o2yT2stWzYm0OQjzuh/h4nHeRbs7p6kF1qGywiIiIirV94wH9nkaeP7kpoA3uYG6q4rIJ7PlzP/LWVj4sFervRJ8yfvmF+9A1rx8DIdnT082jQtY5Myp4JAVsBWUREROQMdut5XflgVSpTXl1O1w4++Hm44ufpWvWrCxf3CaGTf8NCbE37c4r5w1urWJeWw21juhLs686GtBw27M3hxZ0HqXBYnJ0M3905ut51m4+YuzKVBz7dxO3nd+emkdG41rF2dGuggCwiIiJyBgtr58lDV/bho1/SSDlYSG5xGblFZRSUVgAwZ0UqC2eOqHMzk9qsTT3M1DdXUVBSzqwpg7mwd6ej3i8uq2DOij088OlmDhWWNigg7z1cRHGZg398sZX5a/byyLg+DO4ceHIftpkoIIuIiIic4SYNjWTS0KM3GCmvcLBwwz5un7OW91amMiWhc4Ou9fGaNO7+cAMd/dx588ZhtW6E4uHq3KBQfCwnAy9NieP+TzYy/oVlTBoawd1jY2nn5XbS12pKrXduW0REREROmYuzE1f0DyU+OpB/f7WNnMKyesdba3n08y38ce46BkW245NbR5xwl8BTcUGvjnz9P6OYem4M769K4/wnvuejX9JoTQtHKCCLiIiItFHGGP56eS8OF5Xx1OId9Y59aWkSL32fxDXxkbx1YzyB3k03q+vt7sJfLunJpzNGEBHoxVOLd1BS7miy+50stViIiIiItGG9Q/2ZOCSCN5elMDk+kq4dfI4bs2RbBv/4YiuX9gvh71f1abaVJnqF+vHRLcPYn1uMh6tzs9yzITSDLCIiItLG3XlhDzxdnXl44ebj3ks5WMDM99YQ28mPf004+Q1LTpeTk2n05elOlwKyiIiISBvX3sed287vynfbMlmyLaP6eH5JOTe/uQpnJ8OsKYPxcquluaAkD94eD4dTm7HilqWALCIiInIW+P2waKKCvPjbZ5spq3DgcFj+Z+5akg4W8OzkQUQE1rEqxeZPYOc3sOTR07p/K3oG74QUkEVERETOAm4uTvzfpb3YlVnA24m7eebbnXy1+QB/uaQnw7u2P+XrnkG5t8H0kJ6IiIhIW1VeCk7OlT/Ar3p2YETX9vzry20UllYwblAYNwyPOqVLG06+V/lM2GYaNIMsIiIi0vY4KmD5LPhXF/j279WHjTHcd1kvSsod9Av355Ff9z1jQmtz0gyyiIiISFuSvhY+uwPS11S+zk0/6u0enXz5dMYIwgM9W9XSaq2JZpBFRERE2oKSPPj8Hnj5PMjZC+NfhXaRtQ7tFeqHn4drMxd45lBAFhEREWlt9q2HXd81fHzKj/DsUFj+Igy+HmashL4T4BT6hEUtFiIiIiKtR3EufPcwrJgF7r5wz56GnbfkMcDCTd9AeFyTlng2UEAWERERaWnWVq43/MU9/H979x5taV3Xcfz9YVBRUggBUxluhQQVgo2YF3QMMLwUmMalVqagNBaay26wzCSNaGWrhKC4OEiwFCIyGAElTQk1kUFAGSRsnEgOEwyMgNxxZn79cfbInjPnMvuc5znP8+z9fq016+z923s/+7uf3758nt/89m/z0N3wnBfAYw9s+e03rIMdX9TqcFw6tCCcUywkSZKadP//wqeOhH/+Ldh2R3jnv8PPvLnpqrZY6dIvgGwhR5AlSZLqtGED3PiPcMNSWP/DzS+//w7IAnjdKfCyJbBga7j10/Ne5nzoyoxoA7IkSVJd7rkVPvM+GLseXvAS2GHPza+z2yvgVe+H7RfWVERNI7xdSbuzYECWJEmq2pOPjH9x7mtnwjO3hyPOghcfDf4oRycYkCVJkqr0navhyj+AB78HB/wmHPpheNYOzdRiIJ8VA7IkSVJV1n4XPnUU7LQ3vOOz49Mn1DmuYiFJklSVR+4FChx2agXheIC5w0O4kkSTDMiSJEmqXZcyvAFZkiSpbWYzd7ih+cYdyr1brPUBOcmeSZYmubSvbdcky5Kcl+TEJuuTJElqtRYN3XblO4O1BuRegF2TZMWE9sOS3J5k5UwBt5SyqpRy3ITmFwFXllKOBfatuGxJkiTNoCNZd1bqHkE+HzisvyHJAuBM4PWMh9tjkuyb5OeSXDHh385TbPcm4OgkXwS+VGP9kiRJHTbMMbY+tS7zVkq5NsnuE5oPBFaWUlYBJLkYOLyUcirwpi3c9DuAD/W2fynwiYlXSHI8NWz4DAAAEV1JREFUcDzArrvuOrsHIEmSpJHTxBzkFwJ39p0f67VNKslzk5wFHJDkpF7z54D39trvmOx2pZRzSimLSimLdtppp2oqlyRJmk6L5vtq9pr4oZDJxvqnfDaVUtYCSya0rQDeWnFdkiRJFalgasNAYbv9wbz9FT6liYA8BizsO78LsLqBOiRJUpNKgfu+A08+vPllz9wBdthj/muSaCYgLwf2SrIHcBdwNPDrDdQhSZKatPomOPe1U1wY+MOVsO2O81pSe8xmBLr9X8hLB2qEmgNykouAxcCOScYY/2Ld0iQnAFcDC4DzSim31lmHJElqoccfHP97yJ/Bzvs81b7qP+C6M8dHlkc2IFep3skNwzjtuu5VLI6Zov0q4Ko671uSJHXEwpfBbi9/6vyj32+uFm2xdOVXP2ah9b+kJ0mSOuTmi+Chu5uuQhsNcYitkwFZkiRV45H74LIl8Mlfa7oSaU4MyJIkqRob1o3/ffieZutoVJUTcgfY1jBOBG6QAVmSJKlqTm3YTJcyvAFZkiSpbQzYjTIgS5KkhswwpNilIcc2a9N+7EjuNyBLkqRmOVpajcb2Y4sCeEWa+CU9SZJUt3VPwF03Qtmw+WU77AnPef7817SlDMwVcl/OhgFZkqRhdN3fwxdOnvyynfaB372u+vts03/lq3bDHL0NyJIkDaMnHoJsBb952abtX/0Y3PffzdQkdYQBWZKkYZWtYM/XbNr2rX+qLyA7NaLaUfSBtuXofZX8kp4kSW3ywJ3w5b+BJx9tuhLNiQcLE5UOhXhHkCVJaoP16+Dr/wBfOhV++AjsvA/s/fqmq6rXjCOk3QlU1TNgN8mALElS08ZugM+8D+65BX7i5+DuWyZffUKalfYcaHQl9jvFQpKkpjz2AFzxfvj4IfDoWjjyQviVM5quqgFdiU0aFY4gS5JUt0+8Ee5ZsXn7usdh/ZPwsiXwix+AZzwb/u+b819fVVzmrWGTHGjMwxcnh7HbDciSJNXte1+D5+8Huxy4aftWW8N+R8IL9q/+PqdLLbUnmrmGMkeUu2CYFy0xIEuSNB9+8mA4+INNVzFPhnBIUSPFOciSJE3l0e/Dw2uarmIOJhviq3HYb5iHFLdYlQcHA2xrGOc5NMiALEnSRBvWw9fPgY/tBxf/RtPVqIs8WNhchzK8UywkSeq3+ma44n2w+ibY6mnw+ANNVzTEZkhMozwqWmXAbtF+7MpxgyPIkiQBPPEQfO4kOPe18OBd8JalsM+bmqunRaGmdl1JTRoZjiBLknT3LfDJI+Gh/4NFx8LBfwrP3B7+68qmK+uWUQr1neHBx2wYkCVJ+ubF8Mi9cNznYeFLGy5mGALNMDyGDmpoJH4YD4ucYiFJUimw9TY1huMmIsR099nySOOUCzXMgCxJ0nwYqdDX8gCuSmSI/6eg9QE5yRFJzk1yeZLX9dr2TLI0yaVN1ydJatDa78LJ28Hdk/yMsyYP5XVmmpE6CJhClfOwB9qWByVVqjUgJzkvyZokKya0H5bk9iQrk5w43TZKKZeVUt4FvB04qte2qpRyXG2FS5K64bZl439vuaS++/CLZ/UZ6l3rwcJEXeruur+kdz5wBnDBxoYkC4AzgUOBMWB5kmXAAuDUCbc/tpSy8SeM/qR3O0mSpCFXZcBuTzTtyrSMWgNyKeXaJLtPaD4QWFlKWQWQ5GLg8FLKqcBmC04mCfCXwGdLKTfWWa8kSZqDWY+2dyM0dZLTXmaliTnILwTu7Ds/1mubynuAQ4C3JlkCkOS5Sc4CDkhy0mQ3SnJ8khuS3HDvvfdWVLokaThNEewaDRftGfUbnKFM3dbEOsiTvWqmfBcopZwOnD6hbS2wZLo7KaWcA5wDsGjRoi6/y0iS5oMjbZJ6mhhBHgMW9p3fBVjdQB2SJM2PQaYeGNQ1aw39UMgQDkM2EZCXA3sl2SPJ04GjgWUN1CFJqtJ3vwTfvrzpKrTRdKml9kQz1+17kNAFw3wsV/cybxcBXwP2TjKW5LhSyjrgBOBq4DbgklLKrXXWIUmaBxceAZe8rekqalBVmGwiTUx2nzXWMcyJaYtVefAxwLY6MIxbOlDjRnWvYnHMFO1XAVfVed+SJKntZghMHQpUm5nrwUKVBxtd3o8Naf0v6UmSpCHnyPPI6EpXG5AlSXKErRruxxbqSCJtGQOyJEmAQaJK7kt1mwFZkqQpNRj0HI2VGmNAlqRhdNeNsGF901XoRwYJu46+qlu6tDrFljIgS9KwGbsBzn0tXPvXTVeiRk0XWlq+DnJXvsnVNvO834a5lwzIkjRsfnDX+N97bmm2jmFR1ehYE6Fvsvussw6DbbVTYwbaVvtHcbs00GxAliRJzZgxMXUoUW1mrgcLVR5stGc/duUQqtYfCpEkDej6c2HVNZNf9srfg4UHzms5o6M9AaLTZj1E2JXY1EGO6s+KAVmS2uTrZ8PDa2D7hZu233MrbL+bAblO5ogKuTPVbQZkSWqbnzoYfu0Tm7b9xS7N1CJJI8iALElVWfckk/5X/VZPg638yockdYUBWZKqcOOFsOyEyS/bczG87fL5rEZVaXT+pvOipaYYkCWpCt9fBQQO/uCm7bdeBmtXNVKSOqqqUD7dF+ZqX2/LcN+MZg7ohrG3DciSVJWttoaDfn/TtvtWwh1faaYeVaSqj/8mwstk99mmdZBn2LddWjj3R6qseYBtNXFANGB3d6k3DciSRlspcPtV8PiDm1+2zXaw9xtcJmkUdDKItdBs9+MwvsTm+r4xpO876cjjMiBLGm1rboOLf33qy3/nOth5n/mrp1NjLMOmGx/c3eC+bA/7YjYMyJJG27rHx//+8umw52ueal91DXzm92DdE/NfU0dGWCRpWBmQJQngx54HP777U+e33bmxUiRJzXJhTknqDKdfSNJ8MCBLUhc47UKS5o0BWZIqMd3oriO/3VXBgclsV3aY88oaTT4nfc6PkmFcBMaALElVmWyU15HfjpiHNWS3+LlQ4XOm7c/Jqfbtj2rsYPKqMi0OtK0yTd/Wsx8z4HO1S0HagCyp/U7bH856VU0bn+kde57f0bv0CTJs2hQcO2/QfTmM+36uj6mifdKy53W7qpmaq1hIar/7/6f++5j4IdLoh0pXPkIkaTg5gixJkiT1MSBLqs9Dd8P3rmu6CkmSBtL6gJzkiCTnJrk8yev62rdN8o0kb2qyPkmT2LAevn42/N0iOP+N8OSjTVckSdIWqzUgJzkvyZokKya0H5bk9iQrk5w43TZKKZeVUt4FvB04qu+iPwYuqbxoSXOz+mb4+MHw2T+CBVvDhnVQ1jddlSRJW6zuL+mdD5wBXLCxIckC4EzgUGAMWJ5kGbAAOHXC7Y8tpazpnf6T3u1IcgjwbWCbOouXNIAnHoIvngLXnw3P2hHeshR+sBo+/8GmK5ubKlaVqGplikG3M0orYtT6WJvajzXeb+3PjYae80OnZftxiu2ULi7HN4OUmp98SXYHriil/Gzv/MuBk0spv9Q7fxJAKWViON54+wB/CXy+lPKFXtspwLbAvsBjwJtLKRsm3O544Pje2b2B2yt9YKNjR+C+pouQ/dAS9kN72BftYD+0g/0we7uVUnaa2NjEMm8vBO7sOz8GvGya678HOATYLslPlVLOKqV8ACDJ24H7JoZjgFLKOcA5lVU9opLcUEpZ1HQdo85+aAf7oT3si3awH9rBfqheEwF5sgU+pxzGLqWcDpw+xWXnV1STJEmSBDSzisUYsLDv/C7A6gbqkCRJkjbTREBeDuyVZI8kTweOBpY1UIe2jNNU2sF+aAf7oT3si3awH9rBfqhYrV/SS3IRsJjxyeP3AB8qpSxN8gbgY4yvXHFeKeWU2oqQJEmSBlD7KhaSJElSl7T+l/QkSZKk+WRAliRJkvoYkCVJkqQ+BmQNJMmeSZYmuXS6Ns2fJPsmuSTJPyR5a9P1jKokByU5K8nHk/xn0/WMqiSLk3y51xeLm65nlCXZp9cPlyZ5d9P1jCo/o2fHgDxCkpyXZE2SFRPaD0tye5KVSU6cbhullFWllONmatOWqaJPgNcDf1dKeTfwttqKHWIVvTa+XEpZAlwB/GOd9Q6ril4PBXgY2Ibxdfc1CxW9Jm7rvSaOBPyVt1mo63NbM3MVixGS5NWMf3BcUEr52V7bAuA7wKGMf5gsB45hfAm+Uyds4thSypre7S4tpWwyWjlZm6ZXRZ/0/n4IeBR4RSnllfNQ+lCp+LVxCfDOUsoP5qn8oVHR6+G+UsqGJM8D/qaU8hvzVf8wqeo1keRXgBOBM0opn5qv+odF3Z/bmloTPzWthpRSrk2y+4TmA4GVpZRVAEkuBg4vpZwKvGl+Kxw9FfbJ7/beND9dV63DrKp+SLIr8KDheHYqfo+6H3hGHXWOgqr6opSyDFiW5ErAgDwgP7eb4xQLvRC4s+/8WK9tUkmem+Qs4IAkJ03VpjkZtE92T3IOcAHw0ZprGyUD9UPPccAnaqtoNA36evjVJGcDFwJn1FzbqBm0LxYnOb3XH1fVXdwImfPntmbmCLIySduU825KKWuBJTO1aU4G7ZM7gONrq2Z0DdQPAKWUD9VUyygb9PXwafyflLoM2hfXANfUVcwIm/PntmbmCLLGgIV953cBVjdUi8bZJ+1gP7SD/dAe9kU72A/zwICs5cBeSfZI8nTgaGBZwzWNOvukHeyHdrAf2sO+aAf7YR4YkEdIkouArwF7JxlLclwpZR1wAnA1cBtwSSnl1ibrHCX2STvYD+1gP7SHfdEO9kNzXOZNkiRJ6uMIsiRJktTHgCxJkiT1MSBLkiRJfQzIkiRJUh8DsiRJktTHgCxJkiT1MSBLUosl2T3JiqbrmEmSq5Js3zv9cNP1SNJcbN10AZKk7iulvKHpGiSpKo4gS1INklyW5BtJbk1yfF/7w0lOSfLNJNcleV6v/Sd755cn+fBko7BJFiT5aO8630ry2xXUuVOSf+ltc3mSV/baT05yYZIvJvnvJO/qtT8/ybVJbk6yIslBvfY7kuw4Ydvp1bsiyS1Jjuq1L05yTZJLk/xXkk8myVwfiyRVxRFkSarHsaWU7yd5JrA8yb+UUtYC2wLXlVI+kOSvgHcBfw6cBpxWSrkoyZIptnkc8GAp5aVJngF8Ncm/lVL+p/9KSb4MPHuS2/9BKeULE9pOA/62lPKVJLsy/vO1+/Qu2w/4hV7NNyW5EjgGuLqUckqSBcCzptkHvwrsD7wY2LG3H67tXXYA8DPAauCrwCuBr0yzLUmaNwZkSarHe5O8uXd6IbAXsBZ4Erii1/4N4NDe6ZcDR/ROfwr460m2+TpgvyRv7Z3frrfdTQJyKeWgAeo8BNi3bwD3OUk2huvLSymPAY8l+RJwILAcOC/J04DLSik3T7PtVwEXlVLWA/ck+Q/gpcAPgOtLKWMASW4GdseALKklDMiSVLEkixkPni8vpTya5Bpgm97FPyyllN7p9Qz2PhzgPaWUq2e4/0FGkLfq1fnYhG0AlAnXLaWUa5O8GngjcGGSj5ZSLpim3qk80Xd60P0gSbVyDrIkVW874P5eOP5pxqcpzOQ64C2900dPcZ2rgXf3Rm9J8qIk2068UinloFLK/pP8mxiOAf4NOGHjmST79112eJJtkjwXWMz4FIndgDWllHOBpcBLpnlM1wJH9eZO7wS8Grh+mutLUisYkCWpep8Dtk7yLeAjjIffmbwPeH+S64HnAw9Ocp2PA98Gbuwt/XY2cx95fS+wqPelv28D/fOfrweuZLz+j5RSVjMelG9OchPjgf60abb9r8C3gG8CXwT+qJRy9xzrlaTa5an/6ZMkNSXJs4DHSiklydHAMaWUwxus52Tg4VLKZHOhJWmoOedLktrh54EzesudPQAc23A9kjSyHEGWJEmS+jgHWZIkSepjQJYkSZL6GJAlSZKkPgZkSZIkqY8BWZIkSerz/xxWn0HhMNunAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 720x360 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plt.figure(figsize=(10, 5))\n",
"epsilons = np.logspace(-12, 0, 101)\n",
"ax = subplot(111)\n",
"errors = []\n",
"errors2 = []\n",
"axis = np.random.randn(3)\n",
"for epsilon in epsilons:\n",
" #a = np.array([1., 1., 1., np.pi - epsilon])\n",
" a = np.array([axis[0], axis[1], axis[2], epsilon])\n",
" \n",
" a[:3] = a[:3] / np.linalg.norm(a[:3])\n",
"\n",
" R = matrix_from_axis_angle(a)\n",
"\n",
" angle1 = np.arccos((np.trace(R) - 1.0) / 2.0)\n",
"\n",
" r = np.array([R[2, 1] - R[1, 2], R[0, 2] - R[2, 0], R[1, 0] - R[0, 1]])\n",
" angle2 = np.arcsin(np.linalg.norm(r) / 2.0)\n",
"\n",
" error = abs(a[3] - angle1)\n",
" error2 = abs(a[3] - angle2)\n",
" errors.append(error)\n",
" errors2.append(error2)\n",
"plot(epsilons, errors, label=\"via arccos\")\n",
"plot(epsilons, errors2, label=\"via arcsin\")\n",
"\n",
"ax.set_xlabel(\"angle = epsilon\")\n",
"ax.set_ylabel(\"angle error [rad]\")\n",
"ax.set_xscale(\"log\")\n",
"ax.set_yscale(\"log\")\n",
"ax.set_ylim((10e-25, 10e0))\n",
"ax.legend()\n",
"plt.tight_layout()\n",
"plt.savefig(\"small_angles.png\")\n",
"np.mean(errors2), np.median(errors2), np.max(errors2)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.7.4"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment