From 20770146149cfa58f6df19c4c33884accea7bd1d Mon Sep 17 00:00:00 2001 From: Vanessa Vichi Date: Sat, 11 May 2024 10:58:40 +0200 Subject: [PATCH] Conversion from Keplerian elements to attributable elements --- kep_to_att.ipynb | 675 +++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 675 insertions(+) create mode 100644 kep_to_att.ipynb diff --git a/kep_to_att.ipynb b/kep_to_att.ipynb new file mode 100644 index 0000000..557a222 --- /dev/null +++ b/kep_to_att.ipynb @@ -0,0 +1,675 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 1, + "id": "c2e5ae80-979a-454d-b87b-667ad732d63b", + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "from numpy import linalg as LA\n", + "import pandas as pd\n", + "import math\n", + "import scipy.optimize as optimize\n", + "import spiceypy" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "id": "19159865-290b-4ad6-ad47-3a34d79dbb85", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Collecting poliastro\n", + " Using cached poliastro-0.17.0-py3-none-any.whl.metadata (11 kB)\n", + "Requirement already satisfied: astropy<6,>=5.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from poliastro) (5.3.4)\n", + "Collecting astroquery>=0.3.9 (from poliastro)\n", + " Using cached astroquery-0.4.7-py3-none-any.whl.metadata (7.2 kB)\n", + "Requirement already satisfied: jplephem in ./miniconda3/envs/py39/lib/python3.9/site-packages (from poliastro) (2.22)\n", + "Requirement already satisfied: matplotlib!=3.0.1,>=2.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from poliastro) (3.8.0)\n", + "Requirement already satisfied: numba>=0.53.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from poliastro) (0.59.1)\n", + "Requirement already satisfied: numpy in ./miniconda3/envs/py39/lib/python3.9/site-packages (from poliastro) (1.23.5)\n", + "Requirement already satisfied: pandas in ./miniconda3/envs/py39/lib/python3.9/site-packages (from poliastro) (2.2.1)\n", + "Requirement already satisfied: plotly<6,>=4.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from poliastro) (5.22.0)\n", + "Requirement already satisfied: pyerfa in ./miniconda3/envs/py39/lib/python3.9/site-packages (from poliastro) (2.0.1.4)\n", + "Requirement already satisfied: scipy>=1.4.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from poliastro) (1.12.0)\n", + "Requirement already satisfied: PyYAML>=3.13 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from astropy<6,>=5.0->poliastro) (6.0.1)\n", + "Requirement already satisfied: packaging>=19.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from astropy<6,>=5.0->poliastro) (23.2)\n", + "Requirement already satisfied: requests>=2.19 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from astroquery>=0.3.9->poliastro) (2.31.0)\n", + "Requirement already satisfied: beautifulsoup4>=4.8 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from astroquery>=0.3.9->poliastro) (4.12.2)\n", + "Requirement already satisfied: html5lib>=0.999 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from astroquery>=0.3.9->poliastro) (1.1)\n", + "Collecting keyring>=15.0 (from astroquery>=0.3.9->poliastro)\n", + " Using cached keyring-25.2.0-py3-none-any.whl.metadata (20 kB)\n", + "Collecting pyvo>=1.1 (from astroquery>=0.3.9->poliastro)\n", + " Using cached pyvo-1.5.1-py3-none-any.whl.metadata (4.7 kB)\n", + "Requirement already satisfied: contourpy>=1.0.1 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from matplotlib!=3.0.1,>=2.0->poliastro) (1.2.0)\n", + "Requirement already satisfied: cycler>=0.10 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from matplotlib!=3.0.1,>=2.0->poliastro) (0.11.0)\n", + "Requirement already satisfied: fonttools>=4.22.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from matplotlib!=3.0.1,>=2.0->poliastro) (4.25.0)\n", + "Requirement already satisfied: kiwisolver>=1.0.1 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from matplotlib!=3.0.1,>=2.0->poliastro) (1.4.4)\n", + "Requirement already satisfied: pillow>=6.2.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from matplotlib!=3.0.1,>=2.0->poliastro) (10.2.0)\n", + "Requirement already satisfied: pyparsing>=2.3.1 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from matplotlib!=3.0.1,>=2.0->poliastro) (3.0.9)\n", + "Requirement already satisfied: python-dateutil>=2.7 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from matplotlib!=3.0.1,>=2.0->poliastro) (2.8.2)\n", + "Requirement already satisfied: importlib-resources>=3.2.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from matplotlib!=3.0.1,>=2.0->poliastro) (6.1.1)\n", + "Requirement already satisfied: llvmlite<0.43,>=0.42.0dev0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from numba>=0.53.0->poliastro) (0.42.0)\n", + "Requirement already satisfied: tenacity>=6.2.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from plotly<6,>=4.0->poliastro) (8.3.0)\n", + "Requirement already satisfied: pytz>=2020.1 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from pandas->poliastro) (2023.3.post1)\n", + "Requirement already satisfied: tzdata>=2022.7 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from pandas->poliastro) (2023.3)\n", + "Requirement already satisfied: soupsieve>1.2 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from beautifulsoup4>=4.8->astroquery>=0.3.9->poliastro) (2.5)\n", + "Requirement already satisfied: six>=1.9 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from html5lib>=0.999->astroquery>=0.3.9->poliastro) (1.16.0)\n", + "Requirement already satisfied: webencodings in ./miniconda3/envs/py39/lib/python3.9/site-packages (from html5lib>=0.999->astroquery>=0.3.9->poliastro) (0.5.1)\n", + "Requirement already satisfied: zipp>=3.1.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from importlib-resources>=3.2.0->matplotlib!=3.0.1,>=2.0->poliastro) (3.17.0)\n", + "Requirement already satisfied: jaraco.classes in ./miniconda3/envs/py39/lib/python3.9/site-packages (from keyring>=15.0->astroquery>=0.3.9->poliastro) (3.4.0)\n", + "Requirement already satisfied: jaraco.functools in ./miniconda3/envs/py39/lib/python3.9/site-packages (from keyring>=15.0->astroquery>=0.3.9->poliastro) (4.0.1)\n", + "Requirement already satisfied: jaraco.context in ./miniconda3/envs/py39/lib/python3.9/site-packages (from keyring>=15.0->astroquery>=0.3.9->poliastro) (5.3.0)\n", + "Requirement already satisfied: importlib-metadata>=4.11.4 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from keyring>=15.0->astroquery>=0.3.9->poliastro) (7.0.1)\n", + "Requirement already satisfied: SecretStorage>=3.2 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from keyring>=15.0->astroquery>=0.3.9->poliastro) (3.3.3)\n", + "Requirement already satisfied: jeepney>=0.4.2 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from keyring>=15.0->astroquery>=0.3.9->poliastro) (0.8.0)\n", + "Requirement already satisfied: charset-normalizer<4,>=2 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from requests>=2.19->astroquery>=0.3.9->poliastro) (2.0.4)\n", + "Requirement already satisfied: idna<4,>=2.5 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from requests>=2.19->astroquery>=0.3.9->poliastro) (3.4)\n", + "Requirement already satisfied: urllib3<3,>=1.21.1 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from requests>=2.19->astroquery>=0.3.9->poliastro) (1.26.18)\n", + "Requirement already satisfied: certifi>=2017.4.17 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from requests>=2.19->astroquery>=0.3.9->poliastro) (2024.2.2)\n", + "Requirement already satisfied: cryptography>=2.0 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from SecretStorage>=3.2->keyring>=15.0->astroquery>=0.3.9->poliastro) (42.0.5)\n", + "Requirement already satisfied: more-itertools in ./miniconda3/envs/py39/lib/python3.9/site-packages (from jaraco.classes->keyring>=15.0->astroquery>=0.3.9->poliastro) (10.2.0)\n", + "Requirement already satisfied: backports.tarfile in ./miniconda3/envs/py39/lib/python3.9/site-packages (from jaraco.context->keyring>=15.0->astroquery>=0.3.9->poliastro) (1.1.1)\n", + "Requirement already satisfied: cffi>=1.12 in ./miniconda3/envs/py39/lib/python3.9/site-packages (from cryptography>=2.0->SecretStorage>=3.2->keyring>=15.0->astroquery>=0.3.9->poliastro) (1.16.0)\n", + "Requirement already satisfied: pycparser in ./miniconda3/envs/py39/lib/python3.9/site-packages (from cffi>=1.12->cryptography>=2.0->SecretStorage>=3.2->keyring>=15.0->astroquery>=0.3.9->poliastro) (2.21)\n", + "Using cached poliastro-0.17.0-py3-none-any.whl (162 kB)\n", + "Using cached astroquery-0.4.7-py3-none-any.whl (5.3 MB)\n", + "Using cached keyring-25.2.0-py3-none-any.whl (38 kB)\n", + "Using cached pyvo-1.5.1-py3-none-any.whl (910 kB)\n", + "Installing collected packages: pyvo, keyring, astroquery, poliastro\n", + "Successfully installed astroquery-0.4.7 keyring-25.2.0 poliastro-0.17.0 pyvo-1.5.1\n", + "Note: you may need to restart the kernel to use updated packages.\n" + ] + } + ], + "source": [ + "pip install poliastro" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "id": "eb4d9629-fbbe-42cc-bbfe-5fb11aff2302", + "metadata": {}, + "outputs": [], + "source": [ + "from astropy import units as u\n", + "from astropy import constants as const\n", + "from astropy.time import Time\n", + "from astropy.coordinates import get_body_barycentric_posvel\n", + "from poliastro.bodies import Sun, Earth\n", + "from poliastro.twobody import Orbit" + ] + }, + { + "cell_type": "code", + "execution_count": 19, + "id": "c5a3cf10-9e10-4d54-8c08-f8b5cac2b0de", + "metadata": {}, + "outputs": [], + "source": [ + "folder='/home/unipi/v.vichi3/Desktop/'\n", + "df=pd.read_csv(folder+'neos_dataframe.csv')" + ] + }, + { + "cell_type": "code", + "execution_count": 20, + "id": "51bac7c4-0fea-4dd9-b152-6ab13b5ed0b5", + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
Unnamed: 0spkidphaHepoch_mjdeaiomwmamoid
0020000433N10.41604000.22271.45810.83304.28178.90334.730.1500
1120000719N15.59604000.54692.63611.58183.85156.22102.370.2010
2220000887N13.88604000.57102.4729.40110.42350.48289.480.0803
3320001036N9.26604000.53282.66526.69215.50132.48321.690.3450
4420001221N17.38604000.43521.92011.88171.3126.68197.640.1080
\n", + "
" + ], + "text/plain": [ + " Unnamed: 0 spkid pha H epoch_mjd e a i om \\\n", + "0 0 20000433 N 10.41 60400 0.2227 1.458 10.83 304.28 \n", + "1 1 20000719 N 15.59 60400 0.5469 2.636 11.58 183.85 \n", + "2 2 20000887 N 13.88 60400 0.5710 2.472 9.40 110.42 \n", + "3 3 20001036 N 9.26 60400 0.5328 2.665 26.69 215.50 \n", + "4 4 20001221 N 17.38 60400 0.4352 1.920 11.88 171.31 \n", + "\n", + " w ma moid \n", + "0 178.90 334.73 0.1500 \n", + "1 156.22 102.37 0.2010 \n", + "2 350.48 289.48 0.0803 \n", + "3 132.48 321.69 0.3450 \n", + "4 26.68 197.64 0.1080 " + ] + }, + "execution_count": 20, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "df.head()" + ] + }, + { + "cell_type": "code", + "execution_count": 21, + "id": "c68024ad-7917-401b-9f72-3ba3641a7bff", + "metadata": {}, + "outputs": [], + "source": [ + "#Mean anomaly to true anomaly\n", + "def kepler_equation(M,E,e):\n", + " #M=mean anomaly, E=#eccentric anomaly, e=eccentricity\n", + " return E-e*np.sin(E)-M\n", + "def mean_to_true_an(M,e):\n", + " E=optimize.newton(kepler_equation, M, args=(M,e))\n", + " nu=math.atan2(math.sin(E)*math.sqrt(1-e**2),math.cos(E)-e) #true anomaly in radians\n", + " return math.degrees(nu)" + ] + }, + { + "cell_type": "code", + "execution_count": 22, + "id": "f8f1dbef-9764-4732-9cfd-fb5cbebb6e5a", + "metadata": {}, + "outputs": [], + "source": [ + "#Perform the conversion for each element of the DataFrame\n", + "N=df.shape[0]\n", + "nu=np.zeros((N,))\n", + "for el in range(N):\n", + " nu[el]=mean_to_true_an(df['ma'][el],df['e'][el])\n", + "df['nu']=nu" + ] + }, + { + "cell_type": "code", + "execution_count": 23, + "id": "8188f5ba-609c-45f9-89f1-549d9b71fb0c", + "metadata": {}, + "outputs": [], + "source": [ + "#Keplerian to Cartesian elements\n", + "def kep_to_car(a,e,i,om,w,nu):\n", + " orb=Orbit.from_classical(Sun, a<\n", + "\n", + "\n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
xyzvxvyvzEPOCH
0-1.043768-0.948369-0.2671770.007095-0.012608-0.00023760400
10.0098772.220705-0.453875-0.0102860.006398-0.00145060400
2-0.5903150.9149740.038736-0.018393-0.0086000.00335060400
30.8883491.323366-0.282285-0.0082310.011640-0.00716760400
42.6826210.150702-0.1166040.0009930.007884-0.00167160400
\n", + "" + ], + "text/plain": [ + " x y z vx vy vz EPOCH\n", + "0 -1.043768 -0.948369 -0.267177 0.007095 -0.012608 -0.000237 60400\n", + "1 0.009877 2.220705 -0.453875 -0.010286 0.006398 -0.001450 60400\n", + "2 -0.590315 0.914974 0.038736 -0.018393 -0.008600 0.003350 60400\n", + "3 0.888349 1.323366 -0.282285 -0.008231 0.011640 -0.007167 60400\n", + "4 2.682621 0.150702 -0.116604 0.000993 0.007884 -0.001671 60400" + ] + }, + "execution_count": 26, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "car_df.head()" + ] + }, + { + "cell_type": "code", + "execution_count": 27, + "id": "cc6906b5-f2f6-4855-814e-c49c01beb2d8", + "metadata": {}, + "outputs": [], + "source": [ + "#Cartesian to Attributable elements\n", + "def car_to_att(r_sun, v_sun, epoch):\n", + " #r_sun, v_sun are the Cartesian heliocentric elements of the asteroid\n", + " #epoch is in MJD\n", + " earth_posvel=get_body_barycentric_posvel('earth',Time(epoch,format='mjd'))\n", + " r_earth=earth_posvel[0].xyz\n", + " v_earth=earth_posvel[1].xyz\n", + "\n", + " #Cartesian geocentric elements\n", + " r=r_sun-r_earth\n", + " v=v_sun-v_earth\n", + "\n", + " #Rotation to the equatorial reference frame\n", + " roteceq=spiceypy.spiceypy.pxform(\"ECLIPJ2000\",\"J2000\",0)\n", + " r=roteceq.dot(r)\n", + " v=roteceq.dot(v)\n", + "\n", + " rho=LA.norm(r)\n", + " dz=r[0]**2+r[1]**2\n", + " epsilon=np.finfo(float).eps\n", + " if (dz.value<100*epsilon):\n", + " ra=0\n", + " else: ra=np.arctan2(r[1],r[0]) #in radians\n", + " dec=np.arcsin(r[2]/rho)\n", + "\n", + " #Computation of first derivatives\n", + " dadx=np.zeros((3,1))<\n", + "\n", + "\n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
cos_rasin_rasin_decra_dotdec_dotMOID
0-0.084129-0.996455-0.6063470.0066160.0067070.1500
10.3923100.9198330.2296630.0075160.0040550.2010
20.3876450.9218090.4599190.0198810.0116050.0803
30.7910260.6117830.1668710.0112270.0038210.3450
40.9964370.0843450.0237120.0053170.0037980.1080
\n", + "" + ], + "text/plain": [ + " cos_ra sin_ra sin_dec ra_dot dec_dot MOID\n", + "0 -0.084129 -0.996455 -0.606347 0.006616 0.006707 0.1500\n", + "1 0.392310 0.919833 0.229663 0.007516 0.004055 0.2010\n", + "2 0.387645 0.921809 0.459919 0.019881 0.011605 0.0803\n", + "3 0.791026 0.611783 0.166871 0.011227 0.003821 0.3450\n", + "4 0.996437 0.084345 0.023712 0.005317 0.003798 0.1080" + ] + }, + "execution_count": 30, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "attr_df['MOID']=df['moid']\n", + "attr_df.head()" + ] + }, + { + "cell_type": "code", + "execution_count": 31, + "id": "f65e0163-c9d8-492a-9095-0d7d34ffad5b", + "metadata": {}, + "outputs": [], + "source": [ + "attr_df.to_csv(folder+'neos_attr.csv',index=False)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "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.9.19" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}