#!/usr/bin/env python
# coding: utf-8

# # Reading data into Astropy Tables
# 
# ## Objectives
# 
#  - Read ASCII files with a defined format
#  - Learn basic operations with `astropy.tables`
#  - Ingest header information
#  - VOTables

# ## Reading data
# 
# 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


# 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


# In[4]:


# Get the name of the columns
data.colnames


# In[5]:


# Get just the values of a particular column
data['obsid']


# In[6]:


# get the first element
data['obsid', 'redshift'][0]


# Astropy [can read a variety of formats](http://astropy.readthedocs.org/en/stable/io/ascii/index.html#supported-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]


# In[10]:


# the units are also stored, we can extract them too
table['MajDiam'].quantity.to('rad')[0:3]


# In[11]:


# Adding values of different columns
(table['RAh'] + table['RAm'] + table['RAs'])[0:3]


# In[12]:


# adding values of different columns but being aware of column units
(table['RAh'].quantity + table['RAm'].quantity + table['RAs'].quantity)[0:3]


# 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]


# 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]


# 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)


# ## Properties when reading
# 
# the reading of the table has many properties, let's imagine the following easy example:

# 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()


# In[25]:


# We can see the meta as a dictionary, but not as key, value pairs
weather.meta


# 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](http://astropy.readthedocs.org/en/stable/io/ascii/read.html#bad-or-missing-values).
# 

# ## Reading VOTables
# 
# VOTables are an special type of tables which should be self-consistent and can be tied to a particular scheme.
# This mean the file will contain where the data comes from (and which query produced it) and the properties for each field, making it easier to ingest by a machine.

# 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


# In[30]:


# extract one  (NOAA_NUMBER) or all of the columns
NOAA = table.array['NOAA_NUMBER']


# In[31]:


# Show the data
NOAA.data


# In[32]:


# See the mask
NOAA.mask


# In[33]:


# Shee the whole array.
NOAA


# In[34]:


# Convert the table to an astropy table
asttable = table.to_table()


# In[35]:


# See the table
asttable


# 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'))