Our first task with python was to read a csv file using np.loadtxt().
That function has few properties to define the dlimiter of the columns, skip rows, read commented lines, convert values while reading, etc.
However, the result is an array, without the information of the metadata that file may have included (name, units, ...).
Astropy offers a ascii reader that improves many of these steps while provides templates to read common ascii files in astronomy.
In [1]:
from astropy.io import ascii
In [2]:
# Read a sample file: sources.dat
data = ascii.read("sources.dat")
data
Out[2]:
Automatically, read has identified the header and the format of each column. The result is a Table object, and that brings some additional properties.
In [3]:
# Show the info of the data read
data.info
Out[3]:
In [4]:
# Get the name of the columns
data.colnames
Out[4]:
In [5]:
# Get just the values of a particular column
data['obsid']
Out[5]:
In [6]:
# get the first element
data['obsid', 'redshift'][0]
Out[6]:
Astropy can read a variety of formats easily. The following example uses a quite
In [7]:
# Read the data from the source
table = ascii.read("ftp://cdsarc.u-strasbg.fr/pub/cats/VII/253/snrs.dat",
readme="ftp://cdsarc.u-strasbg.fr/pub/cats/VII/253/ReadMe")
In [8]:
# See the stats of the table
table.info('stats')
In [9]:
# If we want to see the first 10 entries
table[0:10]
Out[9]:
In [10]:
# the units are also stored, we can extract them too
table['MajDiam'].quantity.to('rad')[0:3]
Out[10]:
In [11]:
# Adding values of different columns
(table['RAh'] + table['RAm'] + table['RAs'])[0:3]
Out[11]:
In [12]:
# adding values of different columns but being aware of column units
(table['RAh'].quantity + table['RAm'].quantity + table['RAs'].quantity)[0:3]
Out[12]:
In [13]:
# Create a new column in the table
table['RA'] = table['RAh'].quantity + table['RAm'].quantity + table['RAs'].quantity
In [14]:
# Show table's new column
table['RA'][0:3]
Out[14]:
In [15]:
# add a description to the new column
table['RA'].description = table['RAh'].description
In [16]:
# Now it does show the values
table['RA'][0:3]
Out[16]:
In [17]:
# Using numpy to calculate the sin of the RA
import numpy as np
np.sin(table['RA'].quantity)
In [18]:
# Let's change the units...
import astropy.units as u
table['RA'].unit = u.hourangle
In [19]:
# does the sin now works?
np.sin(table['RA'].quantity)
Out[19]:
In [20]:
weather_data = """
# Country = Finland
# City = Helsinki
# Longitud = 24.9375
# Latitud = 60.170833
# Week = 32
# Year = 2015
day, precip, type
Mon,1.5,rain
Tues,,
Wed,1.1,snow
Thur,2.3,rain
Fri,0.2,
Sat,1.1,snow
Sun,5.4,snow
"""
In [21]:
# Read the table
weather = ascii.read(weather_data)
In [22]:
# Blank values are interpreted by default as bad/missing values
weather.info('stats')
In [23]:
# Let's define missing values for the columns we want:
weather['type'].fill_value = 'N/A'
weather['precip'].fill_value = -999
In [24]:
# Use filled to show the value filled.
weather.filled()
Out[24]:
In [25]:
# We can see the meta as a dictionary, but not as key, value pairs
weather.meta
Out[25]:
In [26]:
# To get it the header as a table
header = ascii.read(weather.meta['comments'], delimiter='=',
format='no_header', names=['key', 'val'])
print(header)
When the values are not empty, then the keyword fill_values on read has to be used.
In [27]:
from astropy.io.votable import parse_single_table
In [28]:
# Read the example table from HELIO (hfc_ar.xml)
table = parse_single_table("hfc_ar.xml")
In [29]:
# See the fields of the table
table.fields
Out[29]:
In [30]:
# extract one (NOAA_NUMBER) or all of the columns
NOAA = table.array['NOAA_NUMBER']
In [31]:
# Show the data
NOAA.data
Out[31]:
In [32]:
# See the mask
NOAA.mask
Out[32]:
In [33]:
# Shee the whole array.
NOAA
Out[33]:
In [34]:
# Convert the table to an astropy table
asttable = table.to_table()
In [35]:
# See the table
asttable
Out[35]:
In [36]:
# Different results because quantities are not
print(np.sin(asttable['FEAT_HG_LAT_DEG'][0:5]))
print(np.sin(asttable['FEAT_HG_LAT_DEG'][0:5].quantity))
In [37]:
# And it can also be converted to other units
print(asttable[0:5]['FEAT_AREA_DEG2'].quantity.to('arcmin2'))