From 328b842291e2ca071f2358236c8b2159546aa2aa Mon Sep 17 00:00:00 2001 From: petrush <petrus.hyvonen@sscspace.com> Date: Sat, 21 Mar 2020 22:03:08 +0100 Subject: [PATCH] Intro on Vectors and Coordinates WIP --- examples/Vectors and Coordinates.ipynb | 616 +++++++++++++++++++++++++ 1 file changed, 616 insertions(+) create mode 100644 examples/Vectors and Coordinates.ipynb diff --git a/examples/Vectors and Coordinates.ipynb b/examples/Vectors and Coordinates.ipynb new file mode 100644 index 0000000..8ec7ed7 --- /dev/null +++ b/examples/Vectors and Coordinates.ipynb @@ -0,0 +1,616 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Vectors and Coordinates" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Authors\n", + "Petrus Hyvönen, SSC\n", + "\n", + "## Learning Goals\n", + "* *Define Vectors*: How to define vectors in orekit / hipparchus\n", + "* *Spherical Coordinates*: How to create Spherical Coordinates and transform them between Cartesian and Spherical \n", + "\n", + "## Keywords\n", + "orekit, hipparchus, Vector3D, SphericalCoordinate" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [], + "source": [ + "%matplotlib inline" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Initialize orkit and bring up the python-java interface" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "import orekit\n", + "vm = orekit.initVM()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now set up the pointer to the orekit-data.zip file, using one of the helper files. The file should be in current directory if not specified otherwise." + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "from orekit.pyhelpers import setup_orekit_curdir, absolutedate_to_datetime\n", + "setup_orekit_curdir()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now we are set up to import and use objects from the orekit library. Packages can be imported as they were native Python packages" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [], + "source": [ + "from org.hipparchus.geometry.euclidean.threed import Vector3D, SphericalCoordinates" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [], + "source": [ + "from org.orekit.data import DataProvidersManager, ZipJarCrawler\n", + "from org.orekit.frames import FramesFactory, TopocentricFrame\n", + "from org.orekit.bodies import OneAxisEllipsoid, GeodeticPoint\n", + "from org.orekit.time import TimeScalesFactory, AbsoluteDate, DateComponents, TimeComponents\n", + "from org.orekit.utils import IERSConventions, Constants\n", + "\n", + "from org.orekit.propagation.analytical.tle import TLE, TLEPropagator\n", + "from java.io import File\n", + "\n", + "from math import radians, pi, degrees\n", + "import pandas as pd\n", + "import numpy as np\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Vectors and Coordinates\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Vectors in orekit are based on the `Vector3D` class from the hipparchus library (part of the orekit Python wrapper)." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Note that it needs explicitly to be stated that the input parameters are floats, either by floating point notation or by for example float(10)." + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "<Vector3D: {10; 0; 0}>" + ] + }, + "execution_count": 8, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "a = Vector3D(10.0, 0.0, 0.0)\n", + "a" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "<Vector3D: {-2; 10; 5}>" + ] + }, + "execution_count": 11, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "b = Vector3D(-2.0, 10.0, 5.0)\n", + "b" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "With these vectors a number of fundamental operations can be performed." + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "<Vector3D: {8; 10; 5}>" + ] + }, + "execution_count": 12, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "a.add(b)" + ] + }, + { + "cell_type": "code", + "execution_count": 13, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "<Vector3D: {12; -10; -5}>" + ] + }, + "execution_count": 13, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "a.subtract(b)" + ] + }, + { + "cell_type": "code", + "execution_count": 14, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "<Vector3D: {100; 0; 0}>" + ] + }, + "execution_count": 14, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "a.scalarMultiply(10.0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Note that these are objects that are interfaced through the methods of the class. The normal mathematical + - operators are not available. \n", + "\n", + "_See the documentation for [Hipparcus](https://www.hipparchus.org/apidocs/) and the [Vector3D](https://www.hipparchus.org/apidocs/org/hipparchus/geometry/euclidean/threed/Vector3D.html) class._" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The Vector class has some static helper functions, like the angle between the vectors (in radians):" + ] + }, + { + "cell_type": "code", + "execution_count": 15, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "100.14210615657399" + ] + }, + "execution_count": 15, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "degrees(Vector3D.angle(a,b))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Note that vectors are compared using class comparison methods, not the == sign which compares if the objects are the same!" + ] + }, + { + "cell_type": "code", + "execution_count": 17, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "<Vector3D: {10; 0; 0}>" + ] + }, + "execution_count": 17, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "c = Vector3D(10.0, 0.0, 0.0)\n", + "c" + ] + }, + { + "cell_type": "code", + "execution_count": 18, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "True" + ] + }, + "execution_count": 18, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "c.equals(a)" + ] + }, + { + "cell_type": "code", + "execution_count": 19, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "False" + ] + }, + "execution_count": 19, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "c == a" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Spherical Coordinates " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Quite often in space application we are interested in the spherical coordinates. This is available in the [SphericalCoorindinates](https://www.hipparchus.org/apidocs/org/hipparchus/geometry/euclidean/threed/SphericalCoordinates.html) class." + ] + }, + { + "attachments": { + "image.png": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQQAAADwCAYAAADmfBqxAAAgAElEQVR4Ae19CZhdR3XmkdRabAlLtoxbFggJWVi2ZWwLzDJmkYRkybJlCSOE9z3ElrENeGxDHBaBrXeX1y0RBwgmwASSQALJkExCYEgCWYaBTBaSb8KEkMmXMFlJQiZkAiGExPP9T3Vun6pX211a/bq73ve9796qOnWq7nm3/lf3v6dOEaXPybLAjxHRk5bv/zlZHUjtJAskC4yOBT5BRBuIaLH6fo6IvktE20ani6knyQLJAjNhgXeomcLDM9F4ajNZIFlgdCxwSIHBzxDRgtHpVupJskCywMm2wGYi+kci+t9EtPJkN57aSxZIFhgdCywnot8nom8R0cWj063Uk2SBZIGZsAC/Zbh9JhpPbSYLJAuMjgXuUbzBD49Ol1JPkgWSBWbCAs8jon8hot8homWKSHyKet348ZnoUGozWSBZYOYs8EcWZyR2UDo8c91KLScLJAskCyQLJAskCyQLJAskCyQLzBELFEXxFHzbXM6RI0dWPP7446e10dHv95e31XHkyJFT8zxv5YMBHcePH1/V5lqOHTt2SpZlp7fRceTIkWVtdTz++ONLjx07dkbLfizp9Xqr2+hIdU+uBXDzfl/TJvM8fwa+TeujXlEUT8+yDGspGn+KolibZdnGxgqI6OjRo2eXZXlOSx3jZVk+q42OXq/31KIo4CjW+INBWBTF+Y0VEBEAJc/zLW10AGDzPH92Gx2p7sm1wFOJ6P1Nm0yAoFvu6NGjCRCESRIgCGPMhtM9e/ZsWbly5cezLLunybcoiu/Ht0ldrpPn+SN5nr+J002ORVF8X57nb2lSl+vkef7Goijeyukmx16v11pHnucP53n+tibtc52iKB7KsuztnG5yzPP8wTzPH21Sl+vkef5AlmWPIT05OXnHbBgT87aP/X7/wltvvXV3miGcuAXSI4M+FNIjg26POZ3Ksux5/X7/mUSUHhnUL50AQb/lEyDo9piTqSeffHJBURSXg7hSF5gAQRkiAYJ+yydA0O0x51J4FdXv93cdOXJkibi4BAjKGAkQxF2R3jLoxpgLKfnqCjd7nudXWa4rAYIySgIE/e5IMwTdHrM+lef5i3ERWZZtyrLsuY4LSoCgDJMAQb9DEiDo9pj1KQBCv9//D71e7zzPxSRAUMZJgKDfJQkQdHvM6hTIwzzP7+v3+2cFLiQBgjJQAgT9TkmAoNtj1qbgg16W5b6iKHZGXEQCBGWkaQCEFxLRRyJ+A00kwnUZ7tWIcvU1Ivo2Ef0mEWm/9UlyXUYfvGsdkqei9tOe/ARuBJCHaoYw4BACvUiAoAzUFSAcPnz4IiJ6goh+i4guDdh/qDgACC9SwXA/S0RYM3EmEX1MBbq5hJWdJEDoEdEfK18Wblo7JkDQzHFyE1hUw0QiWpbnnp4kQFDG6QIQLrvssksWL178ZSL6eSI6xWN3Z5EHEPBv/FdE9DdEJFdDYvEQAtvgH3vwOUmAgLYAep8hooWqae2QAEEzx8lL5Hl+sXzNiJYTINSzfweAsGrhwoV/uHjx4j8kolPrtT4l7QGETA38t01JD84Q/g6A8FXOP4mAcJtqew+3LY8JEKQ1TtJ5lmXbsYTYbC4BgmkRf7oDQPgwEf371q1bD/pb8pc6AGGMiP5WDb4LDA0Iow9A+A7nn0RAwGML2v5JblseEyBIa0zDOZYbgx+A6iNHjoyBL3AFBkmAUO8HaAkIL8HAGBsb+5yIh3AfEX1DDZh3EhECheD4l0T0F67eOQABpCEGHuqaH8SSQNk3uaAhIIzDbYWIPoVHgKVLlx5buHAhHk/wBSeygvWLI/LQtnVj4AQIwlLTcQq/go9+9KOLEMkmz/MDR44csT67oe0ECPV+gZaA8AsYGKtXr75XAAI68KAaMNgz89eICLOHq4no1129cwDCo0rPT1nqXa7K/oTLGgACCOi/J6IvEhE4ibsWLVr0uRtvvBEk5gNK/yTrF0e4wQMQEL176JMAYcgk3WYAELBKEceQ5gQIIQvp5QIQbL4btjxWgDLslv0vt9122wYDED6kBsz/JCL8A+ODkHPvVudDBwcgfFLpweBzfSEz+NQEBPTr79RMgF8jfuGss856uYqYhLByaPOvWb84oi7K8Ppz6JMAYcgk3WaUZXnb5OTk1hitMYBw0003bV+2bNmHIdvkm2XZQXyb1OU6WZZdUxTFKznd5FgUBW7eVzWpK+ocyPP8uqVLl3753HPPPcr5a9as+eDKlSs/w2nzuHnz5scwKJYtW/alXq+3HzpYZunSpdgr88nLL7/8Ns4LHXu93r48z2+QcmNjY3i78OTFF1/8RpmP89NPPx1T/CfXrl37AVF2VVEUN4q08/cdHx8f7Nz19Kc//b2Qf+SRR3YsWLDg229605v2F0Vx84MPPrgL+hcsWPBdU9/OnTvvRNmiRYu+bpap9J48z2/Beb/ff0HMfZtkIi0A8rAsy4N4ZIipgh/BJwd9+/btw16OKYSaHlMR/3r4R7+ViN6upve252c273swKIjouBFCbbEi+v6pzm7ajhkC+AG0YYvX+KeqbBt3qOYM4UuqPvwn8IHvxBeFp+K5qhykpvm5SZV9zixAOs0QbFZpmQcAyLLsZTAucwgxKl2AgOXP8GREpOIUIGXKkuKRAZkABZBp2MzGBwaQ/YIaFDcbgIABhkH8G1OthM8cgIApOXSZEa4vVPmYQVR/FDUBAZv8Qjf7TbyaiN4nAGGfKrfxHsdV2btsV5YAwWaVFnkIcZ7n+dV4owA1bQEBZCSm14KMTI5J6vcxACEnov9LRPj3vSvwE/LrwOcagHCzGizvC9TXih2AAMIQg3ZwH4gKEyr/LSKPagLC15UOBgTwG4cFIDAP8pBsQ53DbRr9ut5SlmYINqM0zZuYmFiPmYGs3wYQoK8oCvjXy08CBGUNAQh4TICTz38ioleoxwcfKOD9PwbFWQYg8GB9rTR46NwBCPB8RBvS32Q9EeFxBGChzRxqAgLcoKGb743PE9HzAQgXXXQRHptAmMIe5kwJ6X9Vj0VMRmqXl2YImjmaJxBT3xbPvikglGX5HHgzWnqUAEEZBYDw6KOP4hkdr/jA2L+eiHYQ0dlE9BqL7TiLAWGJAQifVgNNA3VV6VeJ6AOsQB4dgABPQAzax9XgB0H3B8pZqVrDwLqXLFny4zX2ZYAPBQb9bxMR9mGA78SzTj311DctWLDgn9UsCY8m5ueQ6pNzY+AECKbJItLY7Uc9zw+k8zx/qWvDkCaAkOf5ZTZPRtW1BAjKEGKG8ANEhAHHgBD6Ff9MDYxTDEDAasTBzMGiAO/8K0ciWe4ABIjg3/or6h8Zbf4QEa2VddU5dH+rBiCg2vOJ6BPKFwF9/odFixZ9fs2aNX1z9iHaY8BDXesnAYLVLP5M3rlIBUDdOTExgRVs1k8dQJiYmACwXIkfxarsRGYCBGUcAQhsL9jGnCbbTAlnIQyiZxqAYJPlPAxuvPsf+ngAYUjWkXHrggULvl4TEFjVXnaaEhwCl8njZeqarbMcFkyAwJaocQQg9Hq9C8D8P/HEE3hV5fzEAoIiI9/IZKRTYQrDXplGAEKVF3nCLPxNNQABIfCtU+0OAOGZY2Njn2gICJgVYebh28oNb6d+n4h+1zN7GJguAULkHSTF4GgEJxKZ5zqPAQQmI12vHQ3daYagDNICEPC6D6/kPh8JCHgcAZFney7HG4I2ezsOdG/cuPHFDQHhvUR0L0zimCHgTQQITvANeDXr/SRA8JpnuBCuoZgZxG52GgKEiYmJC5iMTIAwbG9fzvnnn7/Lxd346qmyNWD7ly9ffsxwXbZVRfwAdgIaKm8JCAPdNd8yyD4A2LYjwwIIcFhC/APMIKKWdydAkKYNnGNwIwAqcwgB8UGxDxDKstwmd0BOgBBj0Upm3/LlyycMQMA2d74AtVVldbJ6yZIlPx0BCGY9Ld0SEAa6WgACeA3MGm2A8BOKbNX660skQPBZR5Wp8GYvxQ+PrLaA4CIjEyBE/BgnRLBi9NPXXHPNuQYgxL5lqBqKfGSo5G0nMwwIVZcsM4SqLPYkAULAUtg9CY8Ix44dY6+wVoBw5MiRFS4yMgFC4MeYKr4dz80WDiEBQp7DL6HxJwGCx3RlWa7Ba0AOcMKiTWcIuIGxTyPrMY8JEEyLONMgdOFUdHaaIUzZKM0QpmzRyRnWDbCiXq93LhyEOC2PTQABZGRRFE6HEOhPgCCtHD63AAJWD24K15ySSI8MU7bAWZohCHvwgCzL8hIAgijSTusCAtaVS/JQUyYS3L7Isp2m147KKhZAsNnLm5cAQTdPAgRhDwxIMP/Hjh17msgeOo0FBDxqTExMvPb48ePwrQ9+QoCAZdV33333wRQgpQoOc6Asy2tht6ZfM0BKEz22ACkN9EQHSPHo3osAKZ7yGDulACkYqdg9Kc/zh7FOITRyYwABZCSWQff7/e1dBEgBGYkArevWrYMv/HwMkLKOiEr526QZgrSG9bWjLhCRSjMEIlIBUPdikVKEzYJvGZiMhC6fH4LZFpDdzEMaN74gI+frIwN88LVQdBZAwFLjVTYbuvLSI4NumXkPCNh6nQciH3UTDad8MwRzK/e2gNDv9y/Msux5ohfzERAQVXhoUY4FENJrx/TaUQyVmqdZlj233+9XLqltAQHBTEx/9DaAACCwkJHzERCwNwKCi2ifBAiaOWyeirpARGrezhAwUE3ysA0g4HHD1Af7NwEEkJFlWe52bA0/HwFhsNGNeT8nQNAtkvwQdHtEpXy7JzUBBN7K3UVG1gUEJiNxdFzQfAQEqyksgHCYiKw8jFXBCX5mfJavZaguLQFCZQr3iRxYMFhRFPtdrH9dQGAyUgRAHepITUDAPgNXDSnRMxIgKHtYAEG3VEQqkYq6keb8IwMPcixXLooC8facH5Z1CqgCRSpuj5GPBQSQkUVRYBON0Ge+AAICz2ALNecnAYJumjRD0O1hTWHQypgDViGVGTPAIZrn+d6yLK/w6eKyGEBgMjKy/fkCCPcQ0R1sR9sxAYJulQQIuj2sqaIo7rUw9VbZmAGJAV4UxaVdBUgBGVkUxSDwZkz782SjFsRFRIQic08D7XdLgKCZI71l0M2hp8ATFEWxM/afHLV9AxJkJFY+Yit3nx+C3gv3WwYbGelrX+idDzOEt4YeF2APCyDgkSu4qa6wJXQkUlEYZE5yCDLmQOQgG5jEJctkJJOHbQEB0XEALqyPfw9X+1yujvMBEKI8Ri2AkByTkmOSPlzUTYJQWoNP5CBzyuZ5/gw8Jih1g0MbQMArLlefXPmy7XnyyGBcsj2ZAEG3S+IQdHvACch08/U+BhjVh2RBRkKnKdcUELATE2Iymvo4nQCBLRF3TICg2ykBgrAH3Hz7/T5i52ufyEE2qCNlsQzapg+CTQABrzxRT+uckZDtG0UyOZcfGZbKCw2dWwBhIxGtDtWT5YlDkNaYpQFSZEgznGMlIAdA1S/PTxTaZJmMPH78uHPVXB1AmJiYeHG/378aZKTZnpmOAYSnPvWpa5YtW4Z9ALFuovY3z3M4Px1oUlfU2Z9l2TUiXbsfKpT9Qdaxbdu2288888yPcDrmiDgEZVlWOmLqmDK9Xu/KsiwPmfl10iCu8zx/VZ06pmye54hDcJ2ZXzN9eZ7nN9Sso/12WZbtwp4j0IG1PuY9OpJpHjjYXxE31pEjR5a4OsqyrnKZn2XZblcAVCkXCwjwZCzL8uHQ7k6sO9RXrG246667Xol4CHjr0eSLV7D4NqnLdfDKFbEMOd3kCG4GfArXXbBgwc+vWrVqA6djjo899ti6ycnJc2NkXTKPPfbY0yYmJs5zlcfkY6aCx8sYWZcMZip4PHWVx+TjTxGL9WJkXTKK8L4E5S6PXr5fR+aIgaOmi9iSy/sJDTKurG7QqG3AYwCBycgYxyTRB6cfviAj5+IjA/YnxA7OtT6WR4Za9SGcHhl0k83K146YWsVOZ2IAgcnIGFmYLwQIWZZt5a3cuwAE6CqKYrP66eYiIDxGRMFHKv3Wtfoh4G3QEI9k1pPpBAjSGrOQQ8AAy7Isal9FXGpokCMSMpOHIVk2nQ8Qsizbzvog3xYQoA/tcdvpteOUJSwzhOSHMF/8EEAeqsF2VuzAxa3jkoW+fr+/S5KRLtmpW/DEmQ0Q8LyVZdnLTDKyKSDgGQ4rHy1k5FycIZgmjkonQNDNNG9eO7KbL++eFDtwYS6brIuMtMnqJj+RMgFBbeV+tY08bAIIaln1AdOTUfUlAYIyRAIE/e6cF4AA1hP/lPJVY+zAhblMWXUTWclIU1Y391RKAgLIQ6yZmCrVz+oCApORuhYtNScA4eDBg+eNj49/n3ZlNRMWQMDvOuRI5lObOATdOiNNKgpmXet17MBFJSkLfWVZPkdTJhJSVmQPnTIgYBt3rHwcEhAZdQChLMs7mIwUKszTOQEIK1aseOemTZtuNi+uTtoCCHWqD2QTIOgmG1lAMJh1rdexAxeVWBYD0+c2LGW1xiwJAAKceyR5aBEbZMUCQlmWL4Kzj0uPyJ8LgLB2yZIln4xdli6uXTtNgKCZY24tf5ahxi3MunblPMi1TEeiLMuXKDfksxwiVXaMXjy6lGV5vW+mUSmMeMvAZKRCZqcfgtA5FwDhhzdv3nxFAoSpXxWPxmbE7qnSuLM5xSFgMDKzDpLOZ4KYgYv6vBsTk5E+nSgL6QUZCT4DM4MuAqQwGYnrjmlf9X8uAMLTERRmGgBhJRG5gtNaf/70yKCbZWQeGRAfANNwB7Ou9To0cCHMZCRmCFplT8KnF7sx4TUlqjOH4FFVFbkeGSYmJtbjNWUlGAFISnYuAAJsOB2AkPwQ5oIfApj1PM/vkoPDd+4buKiHnZvxTI7zkKxsxyVrkpFtAQFTQxCSsu0afU2AoAxn4RASIMx2QMAzuQqCGvP8PLgVXAMXhdjKXbj5tgYEGxnZBhAQQxELg9Q9rR181yUEZzMgVLsvpRmC+EXVjDZxCMImkYNhUMMlayMjXbKi6epUyiryEDERhsjIJoAAffBXmJiYOLNq0DiR7RtFMjmTgHAfEX2DiJ5cvnz5B1XsAWyz9pdE9Beyk5bzVxDRI5yfAIEtceKYSEXdHq3+ybH8GcuWQfoZahvpNT0jTZ11AQFETcyy6lkACDDFgwoQsPnqrxHRQRUQ9ddNO4n0IiL6JSKqfp9pAgTs5YC2oj+JVNRNNTKkYuRgGPReysLNtyiKl7vISCmrX/pwCrIgD7HvAv7RhyVO5NQBBADBsWPHqhiPLp3Ij+nr7t27L1y2bNmHIdvkm2XZQXyb1EWd008//VMAhFNOOeWrt956637k3XfffbtXr179cZfOc845p3/OOecUshy/GYKKyLwG5wfKsry2Qb3Kdr1eD9dwXUsd+xCYpI0OvL1CYJKWOvYWRXFzSx0I1HILdPT7/Rf47tdpLUMHYhtgWTD1iOziq8eyPhkuy7IM/gUDMpLzbMdYQMjzfEu/378zNtBEqK/gWm644QZEJH6/rV8xeYrAfUaMrEPm9wAIW7duRQyD2A8eNbQ9FqZphhDbn0ouzRAqUwxOZu0MAY5BEW6+Uf+6sATISKC8bh57KgYQQEaC3HS9drRp9gECnKvUe/uZ5BAwJf/OggULvvnoo4/WijtgXm8CBN0iiUPQ7RE9cFGtKIq7MSgNFdakb5BxBQw2bOUeI6vaf7rLMQmPGniTwMuq2wKChYycSUC4CLODsbGxL7qun20aOk4TILyKiLaG2pblaYYgrTF4bF1peyWuS52EVMxgxNRbOTDtie2STy+TkbyVu09WtueaIWCnacUZnMLybQBBbjjD+mY4QAoWI4E/+MiIAkLyQ5jtfgh8o4cGI9x8sZW7cm+uzTdwO3zEFK3f718jychQH7iuDRBARgKsTDKyKSDA6QbRpLlNcZzJGcIEAGHVqlVvjQAEPFI8JPqtnU7TDCEBwnwABNPNN3bg4g60ybrISJusdherhAkI8IzM8/wym2wTQAAZidBuNn0zPEP4NABhw4YN10YAwo8S0RbHNSTXZcMwiUMwDOIajMqLUXPzdckaKgdJUxZkJEJVx8jaZJAnAQFkJADBJVsXELDqM7DoZyZnCF8DIFx55ZWXBADhYiJ6wmUTZcPpWMuAdp/ma9csSxyCbpGRfssgmHWt1+Yg1wqNhJTFvzgGsyFSJaVslWk5YUBgMtIiUmXFAoIiI1/DZGSlYPhkJgFh0Bu+/uGuVTk/FhqY0/TIUHUg9iQBgm6pkQQEC7Ou9Tp24KISZHlZNS5WU2QkYvVi+TOiGzEZaajRkjGAADIyz/OrfaHYhNLZAAjOjXP4OhIgsCVOHNMjg26P6lnfwaxr0rEDF5WwRVZZlrtjnINi9CrPyJsD0/qqvyFAYDJS9TWGLJ0NgFBdv+skAYJumQQIuj0GgOBh1jXpmIGLCoo8vF+r7EmE9GZZtgkyEVPmqhUfIECf3HAm1L5SmgBBGcKy/BkBVs+ujB9xkh4ZdCONzCMDfLA9zLrW65iBw2RkjCwr98lKMrILQFAbcp7PbePoa1/IjSog4DEhZoYzuJRpmiGk145z4bUjmHUsdBE3vfc0NHAkGRmSlQ25ZPEvD09Glm0LCPBklPpYr6t9LlfHUQWE1xHRTUZfnckECLpp0iODsgcIRLUxSfS/i2vggCdAmDO5e5JLVv85TqRMWSYjzd2TmgICL6t2kZFm+7Y+zrAfwqBLlutHDMxfJKKFjj4PZSdA0E2SAEG3R+x0eVDLNnBcZKRN1mi6SkpZRLCFZ6SNjLQMiEqHecIcggK9vdIz0pSV7ZtlIj2KM4SjRBTtTo5rmSZAuJ6InitsFTxNHIJuopHhECIHw6D3pqwil6wxB0xZ/fL1FMvC6aYoih166VSqLiBMTEycx7qntAyfxciM6AzhwPDV+HOmCRD8jVpKEyDoRpn1gMBbueuXNZWKHGSDCpBlMnJKw/BZHUAoiuL248ePXzKsZTgnpq8333zzy2YyQAr6iA1liqJ4Jc6bflOAlCHbpQApckjEDAaWZ9kIN99ajyJFUdwb418QCwh4XJicnLQ+dvC1yCNfl8yT5yAjr7zySrhxz2SAFM11W/avznmaIejWShyCbo9aA1ftxrTbFgDVUBulFzwBvATLsrzCrG9LhwBBkpHMIdj0mHkuQDDIyFHiEKJJRPNapwkQ0B9n6DuzD0inRwbdKrPukQFuvlmWvQFH/VLsKdcgY2lJRoZkuY4PEJiMZPKwLSDgnwPLqlnfqHAI5513HqIov4VtUvc4TYCQ/BDmgh8CbqaYwQjyEMEoY2T5BvXJmmSkT5b14egCBBsZ2QYQHLtfj8QMYenSpb9MROPSLnXOEyDo1kqPDLo9goAg3XxjBy6acMnayEiXrNFVKyCAjIROU7YpICBepGO36hkHhLVr196xatWqSfNa66QTIOjWSoCg28M5cCEGN185OGIHLuraZEFG2rZyt8ka3RwkzRmC9Iw05ZsAAl55og1Tl0rXBQT8i2dEhPDpC0877bR3Lly4EHEN/kbFLFjhaMeZfcopp3zo9ttvv8ApEFGQAEE3UgIE3R7WgQsRgIFJHsYOXNSXsmpZ9eWumANS1uielmRAYDJSekZqghHbwUt5AAseiUzPSClTk0OA9+ffE9EXiQhvJ+5avHjxF6655hoEEnkAwU6IqPY/PV+/0a9ayWkChDOI6NQ6HUmkom6tkSUVDWZd63XswEUllgUJiRiFCKyqKRMJlhVZ1lMMiF6vd0HMbkyxMwTlyfhGQR5a264BCJgZ/J2aCWCg4POFjRs3Xo+9GYgIsSEACH+tyqIPIwwI0dfAggkQ2BInjiMJCBZmXet17MBFJcgq8nCfpsSSiNXb6/VekGXZjRYVQ1kxgIABCrnI9mMfGXI14L9fdWopEX3r9a9/PZZwAxDwTwpA+NehTgcyEiDoBsKMU24wrJfGpdIjg2EnHgwOZl2TZlkt05HANl0y5oBDbJAdoxfEIV4DBmIKVs2EACHLsq284UxM+zVmCF9SA57jR16KRwexcxPiQAIQ/rbqrP9kNRFdA5EECLqhEiDo9ugkhcGAgRGDtJEDhzAYY//JcREhvUxG1hkQPkDAbtWS3Ay1rwwdO0P4lhrwvEfEq4nofQIQMGMCIPg2aZW/7TEiehky6ly/VCDPp4lDwPZytcjO9Mggf5XBGBiNjVqKorgfN5rePXsqNHBAHqrBdlZIVrbgkoU+LKtmMrLOgLABAsjILMteZpKRrvZlH2vMEL5uAMK7ieiwAIQPqXLn3gmi3fVE9DFO17l+rmMepwkQkmPSbHdMwmADs44BZ940rrRv4DAZeezYscE/o0/W1G+TxRbzIA8lGVlnQJiAgA1nEFAV7s0x7ZsyNQDhs2rA82a4nyei5wMQNm3ahGA03yWirxJRzGvHHyGiyseizvVb+j/ISoCgWyZxCLo9gtN1KW4buCiHUQEuABmWd8lyuTyasi4yss6AkICAwYiZgWxTnpvtyzJxHvvI8BI16H9bbZjyDSJ61hlnnPHoggUL8Djxp3KQC/22U21fjDrXb1OGvAQIumUSIOj2aA0ILjIycpANeiNloQ9xFI1uDpJ1BgQDAsjIoihA7Dk/sn2nEFEsIEAFdn/6hPJFAF/wD4sXL/4fp59++mNEhChHjT51rt/VwDQBAvwuznG1actPHIJulZF87ah3cThlDhwfGWnKDmubymFZDGLpGTklceKszoCArsnJSThXBbdP5/bN9ox0HUDgqiDbBuSh4BC4rPaxzvW7lE8TILiac+YnQNBNM+sBAeQhblD9sqZSkYNsUEEtq95mekZOaTtxFjsgFBl5/zve8Y6oRUChvoJ3OHz48KG6AVLWr1//+OrVqz8O/VmWHcQX575vr9d7yT333HOFTSYFSNFt1+v19uV5foPNVjXyUoAUOdBgOJn2nUNWxhwIyfrKuQxkZJ7nDzMZyXBVBW0AACAASURBVPm2YwwggIxUKzNfaovLaNPrswGTkWeeeSb2HagbIOW9RHQv2qwxQwDx+B9t/Yy5fls9mZdmCNIaJ/ivoii0sPy6RDiFpffYKDgs6ZaYlTMEOAbleX4gws03iptgMhIzBLeppkpCAwK7MfFbE+YQpmq7z1yAYOx+3eSRAY8L29FyJCAsJqJfISL2YdA6Hbp+TdiRmCZAgAfm6Y4mrdnpkUE3y6wDBNzQRVHcrV+GO+UaZFwDOzeXZfkipEOyXMc3IEwysi0g4F8jz3PJ8jcBBKxpQL1YQLiPiG7l6zWPvus3ZV3paQKE5Icwn2YIYP1BIMYOXNyMPlls5S49I32y8sZ2DQgbGdkGEBBDsSxLkzXHwIZTUaNP5AzhrUS0yNWA6/pd8rb8BAi6VdJrR90e3oELUenmGztwUc8layMjXbJGV4dcd0EeYumyjYxsAgjQhxiPExMTZ5ptq3/6PySi6yxlwaxIQPDqSYCgmyetZdDt0UnKNRjZzRfPNtyQS5bL5dGUhcchPA9B+kk5nJuyZjmn5YAwPSNZho91AUHGeGQdxhEzhA8Q0Q8S0cNGWTA5xwHhe4lo8PgXNIQSSByCbqmR5hCYWTfdfGMHLi5Vyqqt3F/uIiOlrG4mPcWAAPIwz/O90jNSl6TBAqvYtwzYKYrJSFOPSEsOoUdEeN6P/sxxQIi2AwsmQGBLnDiOLCAYzLrW69iBi0osq7aGZ79+TR8nWJbTriMA4ejRozuYjHTJIT92hoDXRXmeY0Vi6CMBAS7a7yEibGEW9fEAApZD3xajhAExRtYlM00cgqs5Z34CBN00IwkIFmZd63XswEUlyDIZqSmxJGL15nm+x7cmQaqOAQTIwDMysn0JCGgK5N9/JqLY3aGeAVCQfVTnP0FE51nyh7ISIOgmSRyCbo9OUjwYHMy61gbLapmOBF5R4gZ2FGvZMXpBHmJNQhcBUvCogevlZdUx7TtWO+LdO3ZfrngW7cJEwjFDwJqHdwox72kCBN08CRB0e3SSglOQh1nX2ogZOHhuVw5M0bsS+/QyGYmt3OsMCNcMATEeQW5Kz0hf+8IA5gyBi/A49EFOuI4OQPgZIlrjqmPm17l+sy6np+mRAY9cl3EbMcf0yKBbaSQeGcCsw234iSeegIdc8BMaOCAjQdAp9+ZaLtG2xvGOuN/vX8NkZJ0BYQMEkJEAK5OMDF2X6psLEFAMknEQ5sx2HchzAMJpLnlbfp3rt9VH3jQBQnJMmu2OSRgUKvx464GLG80kIyMH2eC+tcnayMg6A8IEBOkZaQ4WW/umjOORgcWwvR12VFrFGebRAQimmDdd5/pdihIg6JZJjkm6Paq3AUa2NekaOLat3F2yNsWmLMjIfr/PQUqrKnUGhAQE0zOyUqhOzPbNcpX2zRAgssO310ICBN2q6ZFBt8dIPDKgS5GDYdB7myxmGbat3G2yugmmUlI2z/PLMPCnSqfOmgAC+nfs2LGnTWkZPpPtD5dWOSFAgOBPIjpSVUOcCEDAbKLR6ro61y+a1k6naYaA2I+8B4XWniuRAEG3zKwHBDxy+MjIyEE2sApkeVk1DKObaipVZ0DgLcLk5CT4jOCOQjF9Xbdu3dpQPIT9+/dft3Llyl+CPvPL8RDWrVv37i1btrzFLI9Jp3gIul1TPISpsdHZGW7EWGUsG+HmW3fmsUctqx4KgCr7FgsIajemB7oiS0FG3n333Qcj4yEgyvLzZL9xjhnCtddei4Cp2Ouxij1pyvnSsdcf0LHWNqPz1THLVMxLc/GXKeZNpxmCbp5ZO0PAzYCt2fTLGU4xeAyX6DmKPLxfz7WnYgYEdqvGq1TJIdi1TeX6+ip2v455ZIDSdUT00SntA6JxHICwYsUKeDfuFGW1TmOuP6Rwmh4ZQs0OlSdA0E0yKwGhKIqbi6KAM03w4xtkXJnJyBhZ1AkNCOwUxWRkF4CADW9FNJ1YQEBX30VEL1e7PCPS8gsACGeddZY1EhLbI3QMXX+oPsqnCRAwy9wY0z7LJEBgS5w4zjpAwO5JeZ4jvFfUJzTIJRkZkuUGfQMCACDJw7aAAA5C6gu8duQuLiSiK7G5KxH9OxF9R+3RcLYgFVm29tF3/bHKpgkQkh/CbPdD4BsoNBiVv8JuxBwIybJOHF2y8GTEykK5e5JLVurDuW1AuMjIpoDAy6otZGRohoBl3diUBWHX/00dcY5NXRcmQNB/zTRD0O0xK2YI7OaLI7ofO3Bdsi4yMlavCQgquOXVtmXOTQBBkZF72TNS/8mi9mXA4qXfJKJvC0D4K2UP1+Imoxl30rx+t6S7JM0QdNskxyTdHs5Bzm6+Ujx24KKOKauY6d1SH5+bspxvHuWAUGQknIGsn7qA4NpwRigPzRBYFG9KsCU8ZgaYIXzZ2NuR5Wof5fXXrqwqTBMg4DEpbfba9Ec5MV5GY7NX22AUzLp2iTZZTUAkpCx2TwIHIYq1UymrFRgJHhAgI6HTKNaSdQChKIo7mYzUlOiJWEDgWngTA0Lx74kIW93N5RkCX3P0MT0y6KYa2UcGg1nXeh07cFGJZQEEoffeLKs1ZkkAEMqyvDakD1VjAQFycPixNGdm1QUE1MeeDINXkAkQdHMmQNDtMZKAYGHWtV7HDlxUUrsxDchITYklEaMXPEGe59dF/JMPWggBAshIrHx8/PHHT4tpP/Itg3l1AAPeDTvNEIR1EiAIY4zaI4OHWdd6HTlwCCRklmVvYDJSU2JJhPSqZdpXHT9+fEMXAVJARqpl2nhNWM1mLF2TWU1mCJVH4hyfIWAdQ9BFXBozAYK0xuAeHBkO4Sr8UzqYda3XoYELYZCHaiu12i7RWmMqIclI5hBscmaea4aAgYkyKR9zXQ1nCFUzcxwQkh/CXPBDALOeZdk91V0bOAkNHElGhmRlUy5Zk4xsCwguMtLVvuxjAoQpayiQlmsZEiDMdkBQqxU3Rw6Gwd3gkzXJSJ/s1K114swmCzLS3Mq9DSDAM9LUx/2wtc9l4hj7yHAWEQ3tVZlmCMKSJ2aS4/hD0nPrpVJMxXr2ipKOHAwDXS5ZkJF4vy0bdMlKGT6XsgqoLucAqCyDYxNAABmJZdrSM1LqxLls3ywT6VhAwEYuQ5uWJEAQlkyAoBtj1EjFod45MsyB4yMjTVmHykE2y4KExEpKBFa1ydcFBLzKQUDV0DJobt/WpsiLAYRnEtFHRJ3qdI4DAkLSVwRqddGek0Qq6sYZydeOeheHU3LgwOXTR0ZK2WFNeg5k1XPpPr1ET9UBhKIoXj4xMRFcpo0WYvq6b9++i0MBUlauXPnpq6666iboM78cIMXMr5NOAVJ0u6YAKfr46CSFGzJWEctGuPlGDTJuF/4FWLrMadcxFhBARvb7/dttaxxsuvm6bGXIy7Js60033bQ9IkDKEHfAOuf4DIEvM/qYZgi6qWbtDAHbwmO3I/1yhlOhQcY1lKfgjZz2HWMAgclI12tHm35fX8Xu1zGPDDb1g7wECLppEiDo9piVgFAUxf0YlPql2FO+QYYaIA/VYIteVu0DBOjDsmomI9sCAmYX2DZOkJEJENRPbXnteB0RBWd48k5JgCCtMVqOScFHBhFzoJPdmJiM5N2TQuDBpnMBAraYB3koycg2gMC7XxtkpAsQ4O1oJUG53zjO8RlC8kOY7X4IfLOGBqOKEXAAnowhWdapBoAVaEBGwpMR/+gsH6vXBgguMrIpIGDg4jUl900cXYBwMxEFd49OgCAsmV476saYLa8dcRNjYHHvYwcu5G2yLjLSJsttyqMJCNCHTV2kDJ83AYQ8z5+NDWVZh3G0AQJmBtixaakhO5RMgKCbJD0y6PYYeQ4BzDoIRNnt2IGLOqYsdBVFsVnq43NTlvPNowQEDHgfuVkXEBA81uXJqPphA4QHiQjPz8HPHAeEZxPR2UEjCIEECMIYoz5DEMy61uvYgYtKUhb6MJg1ZSIhZUX20Cl0YNAqN2S4CDs/sYCgPCNfx2SkUyFZQ6i9A/ESPXWqojkOCNV1xp4kQNAtFTlDwA5kX1GRuP7BcAZ7u8ofhOzTtddIycFoYdY1TVJWK7AkIMtkJGIOWESqrFi9x48f31iW5R1MRlYKLCcxgMBk5NGjR52h2IRq2wxBFPtPEyDo9kmAoNsjEhC40p+pwS9f/2PHM4Ts+10WanTkwehg1jWdLKtlOhLwYlS7MQX/QWP0qmCYtwam9VVvQoAgyciY9tNqx8q0gyXuZVnK1Y5ThZFnCRB0Q9UEhJ9Qg/8WoQVcGgDhNSKv/ikGA/69HMy6pjBy4AxesRVFcbdW2ZMI6eWt3CWH4FE3KPIBgklGhtpXbaUZgjKExQ8BO3Vri9tCv08CBN1CNQEBO51h8GMhHX8KIvqjmFfgXMF6zPP8Fg+zrtWJGThg/UEgxsiycp+s3Mq9C0AAUJhkpK997qOYISC68FUiP+p0jj8yJD+Ek+uHACcwAMJvqJsPe4IgmC92C2v+AbNeluWhWA2hgVOW5YuYPAzJyjZdsiYZ2QYQQB66yEhX+7KPAhB+iohqT5cTIOjWTDME3R41ZwgI9f9Pav8PvPq+j4g+o2usmcIAgVtu5GAYaHfJMhmJi+JuuGS5XB5NWXgcwvMQpJ+UawoIpmek1Ilzs32zXKXxyPDzRHTcUe7NToCgmycBgm6PmoCAygAAzBJeoN48aC4CuvYaqcjBMNBok2UyEm8UZLM2WVkuz6UsPCOxdNkW47EJIGDDmTzP90rPSNk2zmX7ZplIjxPR3xLRFUT0wrrf9evXX4lv3XpSfv369Xs3btyIJeK12+c6GzZsuOJZz3pWKx3r1q3bc955513NOono+4noVpEO9m/dunW7N2/efKBOHVN27dq1u9rqGB8f33n++ecjDH+wzy6Z1atXv+zCCy98has8Jv/MM8/ccdFFFx0UslaHO3E/PqoAAY5x7xX57U4jB8OgEVMWuydhAZCtB6asTYbzWFbtxoQfxvqpCwjYuRmPMVZlIpPbF1nmKWYqHyOiP1cBUB4holrfNWvWvHNsbOzNdetJ+TVr1rxr4cKFb5J5dc870vGDY2Njb6nbtpTvyB5d9ONdHVxL636Mj4//wJIlS44IG+ExwPfBHxNmCN8gIq9Pjk/JUFnEYKjqSFkMNrj6VoXGiZQ1ioaSkGUycqhQZNQBhLIsb0FQVVHdeRroK9jz/0pE8FXALAg/2nvqsLkgMo2dpJ19cRV0oaMoih2woauNmHy5a3eMvE2mi2txOc/Z2nPljUo/eMm+q5+OfDyeY1Phhx3lzbIDg0FTyrKIoRh6D82ymgJHAq8oY27UWEDATTs5Obm/gwApW9VaBdPlGlPMxxyXU2WzY5bkVqrCyBPWEXLu8qnrQgdsGYpN6esDyrgfbezBfJVYlh5qdqgc/eBNeoYKIzO66AceY+WS/cimWQwz3/8Vs46GK0Qd6wxc7MaEm2JiYuLMkPIYvTCqcmCKWlYdAgQmI7GVO9Af+kP9RLmjr3j78jNEtCpGhynj25nalHWlFZ+y38anuOqY+XK1qlkWm3bt2i3qY+bkdUDrwh7MVxnL0kU3wqfmJj3hGsMS3A8Ay3BpXA57ycol+3E1B1J4s/B7RFQtPKxR1y/qGAxDlXBTFEXxUOyPEdILo6rdk8ZCstwZHyAoT8aKjGwBCFiW/QYi+iEiWsxt1zmGdqaO0YW3EriGGFmXTFc6IpzWvH4IPq7J1XczH9cS0Q+zmpY+ifbQ2jUT0kvWLHOkMyMf96aZZ4g0TMYMRjD1mNrEyHI3fLLmDeKTZX04ugDBRkY2BIRlRPQj6r2ubDp0jleSIB1XuzaDCSmQ5VhpCocsmVf33LZata4Oc6McT30nIMAePq7Jo7MqCixLr+R8J138Ll30A5sZuZbse/r/O0S0RfjCYOYaNfv16LQXhQZjnudb4MCE2iFZ2YJL1naDuGSlPpzbAMFFRjYABCzf/RQR2YKjmF2xpS9Zvnz5F7ds2YJXR40/XZBlXeiIWA4ur9EKCF0QkDX7IftUnSuHNITIb/zpoh+4J00v2cgOgSv4rnrLhbc7jR9Vgu35BiMuQMYw8MmaDdlkXTeITdbUh7QJCHmeX4Y8m2wdQLj00kuvXLJkyef37dt3I/pS9zsxMfHSHTt2fODOO++8YeXKlZ/dtm3b9zbRUZbla8uyvKJuXZZHPxD3so0O8ERZlr1ehTiPssX69esff+ELX3i/0Y/Xte1HURRYlr6P9dY9wh5tbap4s1b9UDruzfP8aryds92vI5MHI5udAfuJNwlmjACbrFmX01IW+nxkpJTl+rYjA4IgIyvPSFO+BiDAoQQzA0z7a394Z+oFCxZg2ek/E9GHiGibsVbdqxd8SsyGMj4lTHTFcjw2XS2ILhCvg23vIwhIW9NaXot+VHpGpR8hL9mqw6NyYg5G7J4EJLPFHDBlfdfAsjE/DMv69KEMgKAeYcC8e6dNkYDwOiJ6oil5qMih3arfP6YcRTC1+1ci+gsimiQi7/6FXZFlbQm3BkTX0M9l2GOoPCaD+aoYWZdMFzbtoh+2+KGuPo9MvhyMMAJeA7rcfKVs6AIgixsEW7PFyIZkUI4NZfM8vyFGNgAIiIH4fiL6TSL6BSL6ASL6HiJ6KRHBRTn4MQg36MPg/xcFCvAgwxfpx11h1gwdwTZtAtARu1rVVh95DYkuTV0X12IuS9caiEx0YQ/Jm0U2OySmriXoJTtUcaYzeJBzzAFff1jWJ8NlRVHczGQk57mOMXpBRmZZtivLsg0uPTLfAwjwocAiJfji8wf/4nA2wuvG9xHRzyqZn1NvD+Cp+FH1/fipp57635YtW/Zp5acAGayAxPc7Cgi+TUR/TESbhHcjXhVV4dq7IKm60AE7NSS62HYAlKGduqvCyJMu+jEq9vDFD400x8yJYTDKmAO+nsQMXNTHDZLn+at8umRZSC+TkcwhyLqucwcgXKg8D/EKp/Zn06ZNSzHjMbkVpWijmBUALMywcVir/ont27ePtfBOGzTFHm4xDmKui4QOxb439oNX3NDl4+PjmLU1smlX/ejCpm3tAVubS/Zd9h/Z/KIo7o31sw8NXPy4TEaGZKVBXLJK34vhXQb5loCAxSB4PGg0AEI7UytPvW+hmy6vvZ07d46bG8pIO8Scg3DDvhYNPdwGTXRBdDHXhCMRWV87hq6ni350QUB20Q9wWgAl9Cd03SNZrgbbVXXIKNfAxQXyDcJkpE/WNIhN1kZGtgCE1xLRjxIRbt7anxqEm3Oz1xo6nP0bFaIL1wJQEh2tDQi4ltCydKHfetqVTdv2w7dk39rxUc60DUZXf12yTEbKei5ZKcPnpix+aBsZWRcQ9u7dC7LvnWqFIjdX6wjCLWZnap9SD1mGlZPBdSHQPSpEl8MetQAhhq/y2RNlHpuGqlblXfRDvdFwLtmvGpstJ+Zg9PXbJuu4QRp7NfrY6jqAcMsttzx/bGzsUxs2bLgbKzObfCcmJl6ZZdnuJnW5jk/Hc57znH1Lliz5wkUXXfRKlrcdfTps8o68fQBZR1msfQ6WZTlkj2c/+9mHLr/88m2Rulv3Y1TsgbU4cJzCfekbN7OqzDbIXRdgyuI1oMvzypR16UQ+y4KMzLIM5Jz1UwMQ8Nbgl9euXTt49YkZR51vURRrsyw7CN/1OvWkLHQURfGKkI79+/dvXrJkyc+eccYZr5X1cR6rw6xnpnEtcPE282PTXdgDbcEek5OTW2PbNeW6skfbfqBf8NVhm7Yhd603+kxm8mCM6YOUBXnoIyOlbEg3XDvxz9Pv972EXyQgwO/hk0S0JtSurdxBMCGoKhyPvqYCW8J/wbnmwaHD1hznYYWl9i/TQAfrqo5dEF3cDywprxTXPJHL0mtWrcS5H8xRVQU1TrroR4yXbI0ujZ5onYELWf5hQjdIrF6QkVmWvUGx1V4DRQDC96qBO3Cj9SqzFDo8y+Bc8o9E9FnldYhnfqxshNPR0IpEtQx7J0hbSxNRWV0Qbl0QXbgWOKq1iceg7FEtS48ygCHk+F0MKX8SOvr9/jVtrkUu2fe3NotLYwcuLlEFM9kbY9QYverGvypGFu17AAFLQeFtCJKu0UAESWXpxxlEhL3y/oaIBq8+1U+N0HHwRMSsofp0QVJBx3XXXQcnqWsrxTVPuiC6HPZw9QQrCVebhbZl6aZMKN2FTbvoh7lkP9TvWVtuGQTWa8ENkmXZPdZCS2ZIryQjQ7Ks3gEIcAD6aSKKcmlmXfLo8SxDEAoM/LdJefX6Evlf5fxY5y6Wtx2FDgAcwE3zbrTVMfNcy8FNOV/aYw9XtaG3DOhHv9/Hjk6NP8IebXS07kcofmjjzo1ixZjBiNdu+HFjZPkafbImGemTZX04WgABz/a/RETPk3J1zj2eZVg8hbDrGPhmsFY4nyAfrsqdeKfBS87CySBMOdymo7ZJ8y0Hj7VJw0CsGiCgH5Zrie3CQM7zu0Tr6cIe4Mp8RHd0Z2aLYGgwwgWYf9yQrLxml6yNjHTJSn04NwABy7Z/kYieYcrFpEG4wbnGE7wUpCEG/V9a9GE9Bcq+adtQxiLvzIoguvDWxRvXsQuiK8IezmtgT8WO+mHdpMfXuFnG19ImmCs4IN+SfbPNOZN2DUY2qhwwLlmbMUxZHxlpytr0IU8AAlYmfpiIGjHfkYFHeSMMLFoyP3iT8eTixYu/FsOnmJU5PSpEV6Q9uNu24/esXr16J8fItAnE5HVBhHYRzBVesvizaBNbIuZ6R1LGNhhdkWltsq6LkrK48X1stZR16UP+ww8/vHbFihXwPMxd6wV89VEGwg2znpCcenWJWYDvi9ebjT4tiC4EXxl4N3ZBdNWwh/M61bVg34rGnxb2qNrswh4gum1eslUjc/3EHIxYXoznSNt1m7I2Gc5j2Ri2mmW5ruO4Aq/71q5d+yDkm3zLsrwjy7LbYuqOjY3h7cKTF1988RulPHScffbZn0PZ2rVrPyDLYs+hoyiKW2PlpdyuXbtuP/XUU3/n6quvfht27pZldc/r2MOlu9/v3zRX+oEl+71e79XqWhut3nTcu7MnGxfPvQ1FppWyXMd1hCyTkS4Zzo/Qi9dav0JEzsVDrMt1lDtTu2SM/G+q2UEV8UiQVH+qyvBvXesjdNSqJ4Xvv//+K0455ZT/ouI3yKLo8y4CsfKy9OhGLYJd2KOLfmDWKOOHWro6P7J4MMYYlWVjLJPn+V1MRobkA3ovU+Th+pAeWzmILuw/2YBgQpATPC48xSDLEFMB+ZhBRIfCZk6mQT+qyzKILvhbyCAvlZzvhO3RZucj6MBS37Y6lF+LMy6m7zpQ1kU/1KrfofihobbnbLmK1Lsz5scNDNyBjcSNH7UbEyp59F6vfAzMYCNRvwcHHkWfoiroQn+CgX/48OHTDbJsQgECwmFHfUDMwjUbN3BUBYuQbTm4Rcyb1dIeA92OftxNRNVM09uJE+RwtUlPSNZV3gXxB+9YvCVq4w7t6t+szFdGfTiWTfUM3MH1K5Z4sPVYSFYazCKLfz845ryjKXnYAcGEMGtP3njjjdJjELOUfyIigMVT5DW4zjvoBzxEY4muF7j6obwXnesvXPVkPhb1ANhknjrX/BAs5VVWF/bw9KNqJ3QCm2KGEpKbN+WYKmFRkWUwOm3gkzXZap+s2YAhC6efnySi20y52HQXnmVbtmx5tZoJIEgqBj8G2x8oZ6WhNQy2vtk2pbHJ+fJqxEEAiPaI6F3mBqBd2MO3LJ39EHzXgbIu7BHoR6gLg3Is024b4yKqodkoZAxG7yW4ZG3bhrlkbQ0I2acpz7ztNrmYPDg+4QePkXXJCD7lViL6ivJI/DPlShzlNSh0uJoJ5ivCzdx5OlQPkYzgsDXYqaiLwKOhZekxgNCFPboI5mp6yYaMOe/KxWAMXrtNFj90v98f2ibLJutqQMki6syvEtG5LjlfvkG4+USdZV3oGBGia+OiRYv2dxR4NLgsnYjWubwpu7Cp0uEKbuv8Pc0Cm5esKTPv0w0G7sBmuPHB3rvIyDp6L7jggkeI6ONE1Ihx7ohgWtF29yQH4VbrHjNjU9aqrIS7CDzK/YhZlu7qY5f2aNMPn5esq+/zNr/OwGVZsNWhwcOyAcPiufcNy5Yt+9Cb3/zmbahT94s3ANgDsG49Q/4A9jM08ur25UCe5w+00aGu5SHmdprostnjwQcf3HX77bfX2SPxQFEUDzVpX9TpxKZt+6FebcKZDYuUGi+CC9zHc6cYP2Ds1agf+xmYGYTqROhFBOQPEtFrQrpc5TUIN5cKvPKsdrd2CgUKuiC65HLwQHPOYs+1gFP5DBHBp8P7meZ+eNuWhV30Azoi7kPZbDqvYzC4qMZuGxbQC2IOy3qDwOL6hbrwLIsgy1zNV/ld6OiC6II9Ajsw4U0JFoRhdyrMzIY+LfqBLfCwS9Vgk562y4Vb9KO6plgv2apCOjlhgcDArcwEtrosy0NVRuDEo/ditXvSeQEV1mIQTJj+OXZPstYxM7vSERMH0mzbTIOUjfXoNOsijWtRxK43HqWqCyDABrcgAbUPbIogplpmfOL1Y2Nj8AZt9bugOYBBKLZmqFsAxzY2Demf0+WegTu4btxwzFaHZKWhHLIH1b6JMhyZrOY9Z6KrjWdZVzrAocw00cVkWRt7sI5QjEzfDzM2NvbQjh07HpnpfsAjFZyBXLLv63cqs1jAMXAHkua2YT5ZU7UhOyAPieg9TbdeZ88yAJTZVmyadcTK2+Sgw9i1yCbmzVPxB6JiU7oUqX7sbWMPLEt/4IEHDrWJ6YB+bNq06b1EZF0h6+q/zA8tj5eyrnPXkn2XfMp3WMAYuJWUcg/dV2X41xxIscG50Ivdkz7QZmVeF4E2uyKpECtw6GJrZHRBdHVhD7Es/X4ierfp3RhzScKmWLdyfkwdU0b0wyyKTvuW7EcrSYInLCAGbmUS/Ei2G98mW1UyTpQsAnl8gog01aP7wgAABphJREFUYDFEvckuAm12QVJBR4C0814HCrsgurqwB4KoGtcC70bEpnRukmNeXBc2bRDM1ezGwB0ab3mGClJGMwuYg9zHVpuyvhb37NlziyIPzQClvmpaWVvCDcpAdLUlmEaF6OrCHrgWhz3gbTqp/QCOREc23YGQeI4morIVmTrkJRtVOQnZLcCDPIatZlm7Ji1376JFi37xtttuuxp16n4RsSnP84ezLNtVty7Lj4qOXq+3TTnX7OG+1T0ePXq0tT266EcXNsUgbmsPLNkvy/K1ZVleAVsiUK1296VEcwvAoMw0h1hiyEa0hK3XQTI1+pG6INy60NEl0dUmDgL60XbHIWWPAw3IwyqORMCmeDTEKlXvp0U/Kr1duENXytLJsAXg6prneRRbHQAEAMD72269Hmhj+AKMnC5Iqi50dEF0IXYApvjGJdZKmsvSa1QGGCCI7Isi7BGMh9CiH1WXFdHdOIxepSid2C0AtjrP88P20uFcz2DFNl4/R0TYWKTRpwvCrQsdXRFd2PW5kSFUpS52PrItS6/Zp6ecdtppnxwfH0egGt/rXi8g4Fpg15pta+JduIdrClNCt4DiDGz7GeqCIuUABCxX/mUiavw6rgvPsi504Bm5C6KrretuF9fSRRBV9EPZAwP+JnErmKdOQOiiH13EQTA7nNIOCzgGuVXaIotwWphWjlsrBDI5/mIbz7JR0QGeADv9uJaDB0wxKGZPuzaBWEPL0mP7EdjVylQzBAhd9AN/WtgjoY2butnRlA5YwDLInTUMWWy9/qNE1GjrdfYsa0O4sY4GZFl1jaNCdOFa8jy/uo09OIhqbIzMygjipAubcj8AcEJ1rdMuYjrUajAJn7CAMci9ZlGy+JERZ7Dx1utdEG5d6OiQ6LIFHvXaUhYq8rCx6y904VowQ5F6657XsCneNsENHUvYtU8X/bB5yWqNpMT0WaAOIBw+fHivWpx0TdMedRFoM7ShTEzfRoXo6sIeIDBjl6W7bNPApohYDO/GKn5lF8Fc4Q5t85J19Tvld2yBGoCwadGiRb+6ffv270GdJt+iKO7N8/xVTepynS505Hl+H5Zys84mR2z51VZHR9cyY/04dOjQK5cvX/4bd95551V5nrfuBzb3ybLsRvV7NAqn1/HwmH/qYPyIq8a7X0TyHVpHH1F3sMNOW8INz9Zd6PDFgYy5li6ILlxL252P0A9elh7Tb5tMFzY9dOjQol27dt176aWXbrW1EZOn3nghWG9MTIcYlUmmqQUiAAH7E/x4063Xu/As60LHqBBdXVxLF4TbqPQj1ku26f2d6tW0gAcQsO1YBr6q6e5JXeyw0wVJ1YWOLoiuLuyB+AOYGdT8mTXxLuzRRT/glo3Xm5ghaB1MiZmzgAMQEH/vYwFnFG+nu/As60JHF4SbWPfvvWZfYRfX4lqW7mvXLOuiH55grmZzzrS6lhc5BVLBzFjAAghYEw/2+PlNe9SFZ1kXOw51oQPeekbsgNpm6cIeXfSjK3u03Ta9C/fw2j9CqhBnAQMQgNggDxttvd4F4daVjraEG/qh3G4bE11dXYta9z/j/egiiCps2tY9PO7OTlKNLCAA4Q61e1Kjrde7ILpGRUcXRBcHc+0iEGtoWbrvh+/Cph3ZYwx8AchdX39T2QxboNfr4ZViq63XOyTcGodagxlHheiCPdoGYsW1xC5Ld91CXfwuXfSjC/dw1zWm/G4tsGLp0qU/e+655x7FTKHJ9+jRo9cVRXFnk7pcBzrKsryD002ORVHcDOeYJnVFneuyLLtHpGvbpAt7ZFl2fVmWh9v2o61N0Q8sj2/TjyzLDsLhCDrwqrPb2zdp69ICiGf3aSLa1lQpgne0Jdy6IMu60NEF0QV7wH23qT1Rr4sgql3Yo4t+dOEe3saWqW49C1xPRBvqVdGl8Xyq59RPzSUdbTY7Yct1YY82nAP3Y1SuhfuTjs0sgB2RvowdvojoC0LFbSrvayIvnSYLJAvMEwv8tgIAjrOPpbUACWzKkT7JAskC88wCdykAeISIEBH3T4jozfPMBulykwWSBZQF8K73/xHRl4jo14kIodHTJ1kgWWAeW+CHxSxhHpshXXqyQLIAXn1hdgDeAJumpE+yQLLAPLUAAp58hYguUUc8OiSHkHl6M6TLnt8WeAUR/TkR8SYhb1WzBERITp9kgWSBeWSBG4no34hIRu9FAEw8NvzWPLJDutRkgWQBIvo7Nfi3G9b47yofrx43GWUpmSyQLDDHLPD/ASSl6oCc7XNyAAAAAElFTkSuQmCC" + } + }, + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The conventions used are the mathematical ones, i.e. spherical coordinates are related to Cartesian coordinates as follows:\n", + "\n", + " x = r cos(θ) sin(Φ)\n", + "\n", + " y = r sin(θ) sin(Φ)\n", + "\n", + " z = r cos(Φ)\n", + "\n", + " r = √(x2+y2+z2)\n", + "\n", + " θ = atan2(y, x)\n", + "\n", + " Φ = acos(z/r)\n", + "\n", + "r is the radius, θ (phi) is the azimuthal angle in the x-y plane and Φ (Theta) is the polar (co-latitude) angle. These conventions are different from the conventions used in physics where the meanings of θ and Φ are reversed.\n", + "\n", + "\n", + "\n", + "_image from Wikipedia_" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "There are two ways to construct SphericalCoordinates:" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "From angles `SphericalCoordinates(double r, double theta, double phi)`:\n" + ] + }, + { + "cell_type": "code", + "execution_count": 25, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "<SphericalCoordinates: org.hipparchus.geometry.euclidean.threed.SphericalCoordinates@27adc16e>" + ] + }, + "execution_count": 25, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "s = SphericalCoordinates(10.0, pi, pi/2)\n", + "s" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Note that the printout of SphericalCoordinate is limited. Use a print statement instead, here is a helper function:" + ] + }, + { + "cell_type": "code", + "execution_count": 32, + "metadata": {}, + "outputs": [], + "source": [ + "def print_spherical_coordinate(x):\n", + " print(f\"r={x.r:.2f}, Φ={x.theta:.2f}, θ={x.phi:.2f}\")" + ] + }, + { + "cell_type": "code", + "execution_count": 33, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "r=10.00, Φ=3.14, θ=1.57\n" + ] + } + ], + "source": [ + "print_spherical_coordinate(s)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Or from an cartesian vector `SphericalCoordinates(Vector3D v)`" + ] + }, + { + "cell_type": "code", + "execution_count": 34, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "r=11.36, Φ=1.77, θ=1.11\n" + ] + } + ], + "source": [ + "bs = SphericalCoordinates(b)\n", + "print_spherical_coordinate(bs)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Compare if equal" + ] + }, + { + "cell_type": "code", + "execution_count": 35, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "False" + ] + }, + "execution_count": 35, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "s.equals(bs)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Or convert back to cartesian coordinates:" + ] + }, + { + "cell_type": "code", + "execution_count": 36, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "<Vector3D: {-10; 0; 0}>" + ] + }, + "execution_count": 36, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "s.getCartesian()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Exercises " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The Earth radius is approximately 6378 km at equator. Such constants are available in the orekit module `Constants`. SI units are used." + ] + }, + { + "cell_type": "code", + "execution_count": 43, + "metadata": {}, + "outputs": [], + "source": [ + "from org.orekit.utils import Constants" + ] + }, + { + "cell_type": "code", + "execution_count": 44, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "6378137.0" + ] + }, + "execution_count": 44, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "Constants.WGS84_EARTH_EQUATORIAL_RADIUS" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Assuming a spherical Earth of equatorial radius, what is the cartesian coordinates for an earth centric earth fixed coordinate system for the town of Kiruna, located at 67.8558° N, 20.2253° E? Note the difference in definition of spherical coordinate systems." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Answer:" + ] + }, + { + "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.1" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} -- GitLab