

In [15]:

    
import numpy as np

# Set up matplotlib and use a nicer set of plot parameters
%config InlineBackend.rc = {}
import matplotlib
matplotlib.rc_file("templates/matplotlibrc")
import matplotlib.pyplot as plt
%matplotlib inline

Viewing and manipulating data from FITS tables



In [16]:

    
from astropy.io import fits
##OR 
#import pyfits as fits

Opening the FITS file and viewing table contents

We open the table and list its content.



In [17]:

    
SDSS_tab = fits.open('templates/SDSS_galaxy_tab.fits')



In [18]:

    
SDSS_tab.info()









    



Filename: templates/SDSS_galaxy_tab.fits
No.    Name         Type      Cards   Dimensions   Format
0    PRIMARY     PrimaryHDU      16   (452,)       uint8   
1    /scratch/bernd/sample_data/SDSS_galaxy_data/SDSS_galaxy_tab.fits#1  BinTableHDU     16   868491R x 2C   [E, E]

I'm interested in reading EVENTS, which contains information about each X-ray photon that hit the detector.

To find out what information the table contains, I will print the column names.



In [19]:

    
print(SDSS_tab[1].columns)









    



ColDefs(
    name = 'tot_mass'; format = 'E'
    name = 'sSFR'; format = 'E'
)

Now I'll load the data into a separate variable.



In [20]:

    
SDSS_data = SDSS_tab[1].data

We can extract data from the table by referencing the column name.

For example, I'll make a histogram for the galaxy stellar masses of SDSS.



In [21]:

    
bins = numpy.arange(8,12,0.1)
mass_hist = plt.hist(SDSS_data['tot_mass'], bins=bins)

Making a 2-D histogram with some table data

Next, we create a 2d plot of the specific SFR against the stellar mass for SDSS

Use numpy to make a 2-D histogram and imshow to display it



In [22]:

    
mass = SDSS_data['tot_mass']
sSFR = SDSS_data['sSFR']

select_good = (mass>6.0) & (mass<12.0) & (sSFR>-80)
NBINS = (100,50)

mdensity, yedges, xedges = np.histogram2d(sSFR[select_good],mass[select_good], NBINS, range=[[-13,-9],[9.0,12.0]])

extent = [xedges[0], xedges[-1], yedges[0], yedges[-1]]

plt.imshow(mdensity, extent=extent, interpolation='nearest', origin='lower')

plt.xlabel(r'log(stellar mass/[M${}_\odot$]) ')
plt.ylabel(r'log(sSFR/[M${}_\odot$Gyr${}^{-1}$]')









    Out[22]:





<matplotlib.text.Text at 0x7f41ac65c190>

Close the FITS file



In [23]:

    
SDSS_tab.close()

Creating a FITS table from scratch

First, we create a short arrays that contain data to be stored.



In [24]:

    
objects = ['SDSS J1230+3291','SDSS 1109+3019']
mag = [15.6,numpy.nan]
z = [0.1,0.04]
flag = [0,-100]

Now, we create the columns for the table.



In [25]:

    
columns=[]
columns.append(fits.Column(name='Object',format='15A',unit='',array=objects))
columns.append(fits.Column(name='mag',format='E',unit='mag',array=mag))
columns.append(fits.Column(name='z',format='E',array=z))
columns.append(fits.Column(name='flag',format='I',null=-100,array=flag))

Format definition to be used as follows: FITS format code Description 8-bit bytes L logical (Boolean) 1 X bit * B Unsigned byte 1 I 16-bit integer 2 J 32-bit integer 4 K 64-bit integer 8 A character 1 E single precision floating point 4 D double precision floating point 8 C single precision complex 8 M double precision complex 16 P array descriptor 8 Q array descriptor 16

Next, we prepare and store the FITS file.



In [26]:

    
## before astropy v0.4
table = fits.new_table(columns)
## after astropy v0.4
##table = fits.BinTableHDU.from_columns(columns)

hdulist = fits.HDUList([fits.PrimaryHDU(),table])



In [27]:

    
hdulist.writeto('test_table.fits',clobber=True)









    



WARNING:astropy:Overwriting existing file 'test_table.fits'.






    



WARNING: Overwriting existing file 'test_table.fits'. [astropy.io.fits.file]

Now, we reopen the file to update some information.



In [28]:

    
HDU_update = fits.open('test_table.fits',mode='update')
tab_data=HDU_update[1].data
tab_header = HDU_update[1].header
tab_data['object'][1] = 'SDSS J0012-0010'
tab_data['mag'][1] = 15.0
tab_data['z'][1] = 0.06
tab_data['flag'][1] = 0
## ADD comment entries for each column that can be read  by topcat
tab_header['TCOMM1'] = 'Object names'
tab_header['TCOMM2'] = 'Object magnitudes'
tab_header['TCOMM3'] = 'Object redshift'
tab_header['TCOMM4'] = 'Flag'
HDU_update.flush()



In [28]: