Using NEOS package

With the new poliastro version (0.7.0), a new package is included: NEOs package.

The docstrings of this package states the following:

Functions related to NEOs and different NASA APIs. All of them are coded as part of SOCIS 2017 proposal.

So, first of all, an important question:

What are NEOs?

NEO stands for near-Earth object. The Center for NEO Studies (CNEOS) defines NEOs as comets and asteroids that have been nudged by the gravitational attraction of nearby planets into orbits that allow them to enter the Earth’s neighborhood.

And what does “near” exactly mean? In terms of orbital elements, asteroids and comets can be considered NEOs if their perihelion (orbit point which is nearest to the Sun) is less than 1.3 au = 1.945 * 108 km from the Sun.

In [1]:
%matplotlib inline

from astropy import time
from poliastro.twobody.orbit import Orbit
from poliastro.bodies import Earth
from poliastro.plotting import OrbitPlotter

NeoWS module

This module make requests to NASA NEO Webservice, so you’ll need an internet connection to run the next examples.

The simplest neows function is orbit_from_name(), which return an Orbit object given a name:

In [2]:
from poliastro.neos import neows
In [3]:
eros = neows.orbit_from_name('Eros')

frame = OrbitPlotter()
frame.plot(eros, label='Eros')
Out[3]:
[<matplotlib.lines.Line2D at 0x22c8c1ae278>,
 <matplotlib.lines.Line2D at 0x22c851e01d0>]
../_images/examples_Using_NEOS_package_5_1.png

You can also search by IAU number or SPK-ID (there is a faster neows.orbit_from_spk_id() function in that case, although):

In [4]:
ganymed = neows.orbit_from_name('1036') # Ganymed IAU number
amor = neows.orbit_from_name('2001221') # Amor SPK-ID
eros = neows.orbit_from_spk_id('2000433') # Eros SPK-ID

frame = OrbitPlotter()
frame.plot(ganymed, label='Ganymed')
frame.plot(amor, label='Amor')
frame.plot(eros, label='Eros')
Out[4]:
[<matplotlib.lines.Line2D at 0x22c8c6c60b8>,
 <matplotlib.lines.Line2D at 0x22c8c6c62e8>]
../_images/examples_Using_NEOS_package_7_1.png

Since neows relies on Small-Body Database browser to get the SPK-ID given a body name, you can use the wildcards from that browser: * and ?.

Keep it in mind that orbit_from_name() can only return one Orbit, so if several objects are found with that name, it will raise an error with the different bodies.
In [5]:
neows.orbit_from_name('*alley')
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-5-59e3c2e6ee39> in <module>()
----> 1 neows.orbit_from_name('*alley')

c:\users\antonio\desktop\proyectos\poliastro\src\poliastro\neos\neows.py in orbit_from_name(name, api_key)
    126
    127     """
--> 128     spk_id = spk_id_from_name(name)
    129     if spk_id is not None:
    130         return orbit_from_spk_id(spk_id, api_key)

c:\users\antonio\desktop\proyectos\poliastro\src\poliastro\neos\neows.py in spk_id_from_name(name)
    101         for body in object_list[:obj_num]:
    102             bodies += body.string + '\n'
--> 103         raise ValueError(str(len(object_list)) + ' different bodies found:\n' + bodies)
    104     else:
    105         raise ValueError('Object could not be found. You can visit: ' +

ValueError: 6 different bodies found:
903 Nealley (1918 EM)
2688 Halley (1982 HG1)
14182 Alley (1998 WG12)

Note that epoch is provided by the Web Service itself, so if you need orbit on another epoch, you have to propagate it:
In [6]:
eros.epoch.iso
Out[6]:
'2017-09-04 00:00'
In [7]:
epoch = time.Time(2458000.0, scale='tdb', format='jd')
eros_november = eros.propagate(epoch - eros.epoch)
eros_november.epoch.iso
Out[7]:
'2017-09-03 12:00'

Given that we are using NASA APIs, there is a maximum number of requests. If you want to make many requests, it is recommended getting a NASA API key. You can use your API key adding the api_key parameter to the function:

In [8]:
neows.orbit_from_name('Toutatis', api_key='DEMO_KEY')
Out[8]:
1 x 4 AU x 0.4 deg orbit around Sun (☉)

DASTCOM5 module

This module can also be used to get NEOs orbit, in the same way that neows, but it have some advantages (and some disadvantages).

It relies on DASTCOM5 database, a NASA/JPL maintained asteroid and comet database. This database has to be downloaded at least once in order to use this module. According to its README, it is updated typically a couple times per day, but potentially as frequently as once per hour, so you can download it whenever you want the more recently discovered bodies. This also means that, after downloading the file, you can use the database offline.

The file is a ~230 MB zip that you can manually download and unzip in ~/.poliastro or, more easily, you can use

dastcom5.download_dastcom5()

The main DASTCOM5 advantage over NeoWs is that you can use it to search not only NEOs, but any asteroid or comet. The easiest function is orbit_from_name():

In [9]:
from poliastro.neos import dastcom5
In [10]:
atira = dastcom5.orbit_from_name('atira')[0] # NEO
wikipedia = dastcom5.orbit_from_name('wikipedia')[0] # Asteroid, but not NEO.
frame = OrbitPlotter()
frame.plot(atira, label='Atira (NEO)')
frame.plot(wikipedia, label='Wikipedia (asteroid)')
Out[10]:
[<matplotlib.lines.Line2D at 0x22c8cdf2828>,
 <matplotlib.lines.Line2D at 0x22c8cdf2a58>]
../_images/examples_Using_NEOS_package_19_1.png

Keep in mind that this function returns a list of orbits matching your string. This is made on purpose given that there are comets which have several records in the database (one for each orbit determination in history) what allow plots like this one:

In [11]:
halleys = dastcom5.orbit_from_name('1P')

frame = OrbitPlotter()
frame.plot(halleys[0], label='Halley')
frame.plot(halleys[5], label='Halley')
frame.plot(halleys[10], label='Halley')
frame.plot(halleys[20], label='Halley')
frame.plot(halleys[-1], label='Halley')
Out[11]:
[<matplotlib.lines.Line2D at 0x22c8cb79668>,
 <matplotlib.lines.Line2D at 0x22c8cb34e80>]
../_images/examples_Using_NEOS_package_21_1.png

While neows can only be used to get Orbit objects, dastcom5 can also provide asteroid and comet complete database. Once you have this, you can get specific data about one or more bodies. The complete databases are ndarrays, so if you want to know the entire list of available parameters, you can look at the dtype, and they are also explained in documentation API Reference:

In [12]:
ast_db = dastcom5.asteroid_db()
comet_db = dastcom5.comet_db()
ast_db.dtype.names[:20] # They are more than 100, but that would be too much lines in this notebook :P
Out[12]:
('NO',
 'NOBS',
 'OBSFRST',
 'OBSLAST',
 'EPOCH',
 'CALEPO',
 'MA',
 'W',
 'OM',
 'IN',
 'EC',
 'A',
 'QR',
 'TP',
 'TPCAL',
 'TPFRAC',
 'SOLDAT',
 'SRC1',
 'SRC2',
 'SRC3')
Asteroid and comet parameters are not exactly the same (although they are very close):

With these ndarrays you can classify asteroids and comets, sort them, get all their parameters, and whatever comes to your mind.

For example, NEOs can be grouped in several ways. One of the NEOs group is called Atiras, and is formed by NEOs whose orbits are contained entirely with the orbit of the Earth. They are a really little group, and we can try to plot all of these NEOs using asteroid_db():

Talking in orbital terms, Atiras have an aphelion distance, Q < 0.983 au and a semi-major axis, a < 1.0 au. Visiting documentation API Reference, you can see that DASTCOM5 provides semi-major axis, but doesn’t provide aphelion distance. You can get aphelion distance easily knowing perihelion distance (q, QR in DASTCOM5) and semi-major axis Q = 2*a - q, but there are probably many other ways.

In [13]:
aphelion_condition = 2 * ast_db['A'] - ast_db['QR'] < 0.983
axis_condition = ast_db['A'] < 1.3
atiras = ast_db[aphelion_condition & axis_condition]

The number of Atira NEOs we use using this method is:

In [14]:
len(atiras)
Out[14]:
16

Which is consistent with the stats published by CNEOS

Now we’re gonna plot all of their orbits, with corresponding labels, just because we love plots :)

In [15]:
from poliastro.twobody.orbit import Orbit
from poliastro.bodies import Earth

earth = Orbit.from_body_ephem(Earth)

We only need to get the 16 orbits from these 16 ndarrays.

There are two ways:

  • Gather all their orbital elements manually and use the Orbit.from_classical() function.
  • Use the NO property (logical record number in DASTCOM5 database) and the dastcom5.orbit_from_record() function.

The second one seems easier and it is related to the current notebook, so we are going to use that one:

We are going to use ASTNAM property of DASTCOM5 database:

In [16]:
frame = OrbitPlotter()

frame.plot(earth, label='Earth')

for record in atiras['NO']:
    ss = dastcom5.orbit_from_record(record)
    frame.plot(ss)
../_images/examples_Using_NEOS_package_35_0.png
This slightly incorrect, given that Earth coordinates are in a different frame from asteroids. However, for the purpose of this notebook, the effect is barely noticeable.

If we needed also the names of each asteroid, we could do:

In [17]:
frame = OrbitPlotter()

frame.plot(earth, label='Earth')

for i in range(len(atiras)):
    record = atiras['NO'][i]
    label = atiras['ASTNAM'][i].decode().strip() # DASTCOM5 strings are binary
    ss = dastcom5.orbit_from_record(record)
    frame.plot(ss, label=label)
../_images/examples_Using_NEOS_package_38_0.png
We knew beforehand that there are no Atira comets, only asteroids (comet orbits are usually more eccentric), but we could use the same method with com_db if we wanted.

Finally, another interesting function in dastcom5 is entire_db(), which is really similar to ast_db and com_db, but it returns a Pandas dataframe instead of a numpy ndarray. The dataframe has asteroids and comets in it, but in order to achieve that (and a more manageable dataframe), a lot of parameters were removed, and others were renamed:

In [18]:
db = dastcom5.entire_db()
db.columns
Out[18]:
Index(['NUMBER', 'NOBS', 'OBSFRST', 'OBSLAST', 'EPOCH', 'CALEPOCH', 'MA', 'W',
       'OM', 'IN', 'EC', 'A', 'QR', 'TP', 'TPCAL', 'TPFRAC', 'SOLDAT', 'DESIG',
       'IREF', 'NAME'],
      dtype='object')

Also, in this function, DASTCOM5 data (specially strings) is ready to use (decoded and improved strings, etc):

In [19]:
db[db.NAME == 'Halley'] # As you can see, Halley is the name of an asteroid too, did you know that?
Out[19]:
NUMBER NOBS OBSFRST OBSLAST EPOCH CALEPOCH MA W OM IN EC A QR TP TPCAL TPFRAC SOLDAT DESIG IREF NAME
2687 2688 1732 19500814 20160827 2455701.5 20110520.0 22.786977 183.484408 95.422520 3.454589 0.143215 3.165183 2.711882 2.455571e+06 2.011011e+07 0.309024 2.457850e+06 1982 HG1 23 Halley
737341 900001 161 0 0 1633920.5 -2390607.0 0.165294 88.110000 30.810000 163.470000 0.967600 18.067901 0.585400 1.633908e+06 -2.390525e+06 0.620000 0.000000e+00 1P SAO/-239 Halley
737342 900002 161 0 0 1661840.5 -1631115.0 0.031113 89.110000 32.060000 163.700000 0.967700 18.095975 0.584500 1.661838e+06 -1.631113e+06 0.070000 0.000000e+00 1P SAO/-163 Halley
737343 900003 161 0 0 1689880.5 -860823.0 0.211345 90.778000 34.018000 163.340000 0.967680 18.118812 0.585600 1.689864e+06 -8.608065e+05 0.962000 0.000000e+00 1P SAO/-86 Halley
737344 900004 161 0 0 1717320.5 -111008.0 359.963217 92.559000 35.904000 163.589000 0.967370 17.995709 0.587200 1.717323e+06 -1.110108e+05 0.349000 0.000000e+00 1P SAO/-11 Halley
737345 900005 161 0 0 1745200.5 660206.0 0.142118 92.652000 36.129000 163.577000 0.967550 18.030817 0.585100 1.745189e+06 6.601260e+05 0.460000 0.000000e+00 1P SAO/66 Halley
737346 900006 161 0 0 1772640.5 1410324.0 0.019990 93.694000 37.219000 163.437000 0.967840 18.132463 0.583140 1.772639e+06 1.410322e+06 0.934000 0.000000e+00 1P SAO/141 Halley
737347 900007 161 0 0 1800800.5 2180429.0 359.761537 94.147000 37.908000 163.574000 0.967980 18.159588 0.581470 1.800819e+06 2.180518e+06 0.223000 0.000000e+00 1P SAO/218 Halley
737348 900008 161 0 0 1828920.5 2950425.0 0.057332 95.241000 39.111000 163.367000 0.968750 18.429120 0.575910 1.828916e+06 2.950420e+06 0.898000 0.000000e+00 1P SAO/295 Halley
737349 900009 161 0 0 1857720.5 3740301.0 0.158377 96.510000 40.579000 163.542000 0.968590 18.375995 0.577190 1.857708e+06 3.740216e+06 0.842000 0.000000e+00 1P SAO/374 Halley
737350 900010 161 0 0 1885960.5 4510625.0 359.959606 97.028000 41.210000 163.479000 0.968910 18.454165 0.573740 1.885964e+06 4.510628e+06 0.749000 0.000000e+00 1P SAO/451 Halley
737351 900011 161 0 0 1914920.5 5301008.0 0.135791 97.582000 41.974000 163.394000 0.968710 18.395334 0.575590 1.914910e+06 5.300927e+06 0.630000 0.000000e+00 1P SAO/530 Halley
737352 900012 161 0 0 1942840.5 6070318.0 0.032109 98.799000 43.261000 163.476000 0.968040 18.173655 0.580830 1.942838e+06 6.070315e+06 0.976000 0.000000e+00 1P SAO/607 Halley
737353 900013 161 0 0 1971160.5 6840929.0 359.952171 99.149000 43.800000 163.418000 0.968150 18.197174 0.579580 1.971164e+06 6.841003e+06 0.267000 0.000000e+00 1P SAO/684 Halley
737354 900014 161 0 0 1998800.5 7600602.0 0.157833 99.997000 44.687000 163.443000 0.967850 18.097667 0.581840 1.998788e+06 7.600521e+06 0.171000 0.000000e+00 1P SAO/760 Halley
737355 900015 161 0 0 2026840.5 8370310.0 0.124640 100.101000 44.930000 163.447000 0.967810 18.090090 0.582320 2.026831e+06 8.370228e+06 0.770000 0.000000e+00 1P SAO/837 Halley
737356 900016 161 0 0 2054360.5 9120714.0 359.940520 100.777000 45.646000 163.311000 0.968070 18.169746 0.580160 2.054365e+06 9.120719e+06 0.174000 0.000000e+00 1P SAO/912 Halley
737357 900017 161 0 0 2082520.5 9890819.0 359.774025 101.484000 46.561000 163.399000 0.967890 18.122392 0.581910 2.082538e+06 9.890906e+06 0.188000 0.000000e+00 1P SAO/989 Halley
737358 900018 161 0 0 2110480.5 10660308.0 359.839206 102.473000 47.624000 163.112000 0.968870 18.454867 0.574500 2.110493e+06 1.066032e+07 0.434000 0.000000e+00 1P SAO/1066 Halley
737359 900019 161 0 0 2139360.5 11450402.0 359.793477 103.704000 49.054000 163.224000 0.968790 18.416854 0.574790 2.139377e+06 1.145042e+07 0.061000 0.000000e+00 1P SAO/1145 Halley
737360 900020 161 0 0 2167680.5 12221015.0 0.201554 103.849000 49.304000 163.192000 0.968840 18.427792 0.574210 2.167664e+06 1.222093e+07 0.323000 0.000000e+00 1P SAO/1222 Halley
737361 900021 161 0 0 2196560.5 13011109.0 0.179564 104.500000 50.152000 163.076000 0.968930 18.432893 0.572710 2.196546e+06 1.301103e+07 0.082000 0.000000e+00 1P SAO/1301 Halley
737362 900022 161 0 0 2224680.5 13781105.0 359.927910 105.295000 51.020000 163.113000 0.968370 18.216883 0.576200 2.224686e+06 1.378111e+07 0.187000 0.000000e+00 1P SAO/1378 Halley
737363 900023 161 0 0 2253040.5 14560628.0 0.234781 105.835000 51.866000 162.890000 0.968000 18.115625 0.579700 2.253022e+06 1.456061e+07 0.133000 0.000000e+00 1P SAO/1456 Halley
737364 900024 161 0 0 2280480.5 15310814.0 359.842327 106.976000 53.057000 162.917000 0.967750 18.021705 0.581200 2.280493e+06 1.531083e+07 0.739000 0.000000e+00 1P SAO/1531 Halley
737365 900025 161 0 0 2308300.5 16071024.0 359.954121 107.550300 53.770100 162.905500 0.967490 17.951861 0.583615 2.308304e+06 1.607103e+07 0.040600 0.000000e+00 1P H593/0 Halley
737366 900026 161 0 0 2308300.5 16071024.0 359.954121 107.550300 53.770100 162.905500 0.967490 17.951861 0.583615 2.308304e+06 1.607103e+07 0.040600 0.000000e+00 1P SAO/1607 Halley
737367 900027 278 0 0 2335640.5 16820831.0 359.805545 109.221400 55.566900 162.264900 0.967933 18.168865 0.582621 2.335656e+06 1.682092e+07 0.779400 0.000000e+00 1P I353/0 Halley
737368 900028 278 0 0 2335640.5 16820831.0 359.805545 109.221400 55.566900 162.264900 0.967933 18.168865 0.582621 2.335656e+06 1.682092e+07 0.779400 0.000000e+00 1P SAO/1682 Halley
737369 900029 718 0 0 2363600.5 17590321.0 0.101706 110.709300 57.245800 162.372500 0.967686 18.087083 0.584466 2.363593e+06 1.759031e+07 0.562300 0.000000e+00 1P J103/0 Halley
737370 900030 718 0 0 2363600.5 17590321.0 0.101706 110.709300 57.245800 162.372500 0.967686 18.087083 0.584466 2.363593e+06 1.759031e+07 0.562300 0.000000e+00 1P SAO/1759 Halley
737371 900031 653 0 0 2391600.5 18351118.0 0.020156 110.704300 57.518500 162.258800 0.967394 17.989419 0.586563 2.391599e+06 1.835112e+07 0.939600 0.000000e+00 1P SAO/1835 Halley
737372 900032 653 0 0 2418800.5 19100509.0 0.243754 111.737100 58.562900 162.218600 0.967302 17.958530 0.587208 2.418782e+06 1.910042e+07 0.678500 0.000000e+00 1P SAO/1910 Halley
737373 900033 7428 18350821 19940111 2449400.5 19940217.0 38.384264 111.332485 58.420081 162.262691 0.967143 17.834144 0.585978 2.446467e+06 1.986021e+07 0.395317 2.452124e+06 1P J863/77 Halley

Panda offers many functionalities, and can also be used in the same way as the ast_db and comet_db functions:

In [20]:
aphelion_condition = (2 * db['A'] - db['QR']) < 0.983
axis_condition = db['A'] < 1.3
atiras = db[aphelion_condition & axis_condition]
In [21]:
len(atiras)
Out[21]:
347

What? I said they can be used in the same way!

Dont worry :) If you want to know what’s happening here, the only difference is that we are now working with comets too, and some comets have a negative semi-major axis!

In [22]:
len(atiras[atiras.A < 0])
Out[22]:
331

So, rewriting our condition:

In [23]:
axis_condition = (db['A'] < 1.3) & (db['A'] > 0)
atiras = db[aphelion_condition & axis_condition]
len(atiras)
Out[23]:
16