In [18]:

    

import numpy as np
from matplotlib import pyplot as plt
from astropy.io import fits
from GetR200 import getr200
%matplotlib inline

altnames={'NRGb004':'USGU189', 'NRGs027':'HDCE 0511', 'NRGs038':'', 'NRGs076':'', 'NRGs090':'USGCU309', 'NRGs110':'Abell1137',
       'NRGs117':'', 'NRGb128':'HCG 051', 'NRGb155':'Abell1367', 'NRGb177':'', 'NRGb226':'Coma', 'NRGb244':'',
       'NRGb247':'MKW11', 'NRGs317':'', 'Abell2063':''}


    

In [31]:

    

class cluster:
    def __init__(self,name,vr,sigma,ra,dec,r200):
        self.name=name
        self.vr=vr
        self.sigma=sigma
        self.ra=ra
        self.dec=dec
        self.r200=r200
    def readNSA(self):
        infile=self.name+'_NSA.fits'
        self.nsa=fits.getdata(infile)
        #infile.close()
    def readAGC(self):
        infile=self.name+'_AGC.fits'
        self.agc=fits.getdata(infile)
        #infile.close()
    def plotphasespace(self):
        deltaR = np.sqrt((self.ra - self.agc.RA)**2+(self.dec-self.agc.DEC)**2)
        deltaV = (self.agc.VOPT - self.vr)/self.sigma
        HIflag = self.agc.FLUX100 > 0
        plt.figure(figsize=(6,4))
        plt.plot(deltaR/self.r200,deltaV,'ko',label='AGC')
        plt.plot(deltaR[HIflag]/self.r200,deltaV[HIflag],'bs',markersize=10,mec='b',mfc='None',label='HI')
        deltaR = np.sqrt((self.ra - self.nsa.RA)**2+(self.dec-self.nsa.DEC)**2)
        deltaV = (self.nsa.ZDIST*3.e5 - self.vr)/self.sigma
        plt.plot(deltaR/self.r200,deltaV,'r.',label='NSA')
        plt.axis([0,3,-3,3])
        plt.title(self.name+' ('+altnames[self.name]+')',fontsize=20)
        plt.xlabel('$\Delta r/R_{200} $',fontsize=20)
        plt.ylabel('$\Delta v/\sigma $',fontsize=20)
        plt.legend(loc='upper right',numpoints=1)
        plt.savefig(self.name+'_phasespace.png')


    

In [23]:

    

##### BEGINNING OF MAIN PROGRAM  ######

# read in sample.dat
infile=open('sample.dat','r')
ra=[]
dec=[]
vr=[]
name=[]
for line in infile:
    #print line
    t=line.split()
    name.append(t[0])
    ra.append(float(t[1]))
    dec.append(float(t[2]))
    vr.append(float(t[3]))
infile.close()
# convert the lists into an array
ra=np.array(ra,'f')
dec=np.array(dec,'f')
vr=np.array(vr,'f')
name=np.array(name,'S10')


    

In [24]:

    

# read in velocity dispersions from biweight_center_scale.dat
infile= open('biweight_center_scale.dat','r')
name2=[]
centerv2=[]
sigma2=[]
for line in infile:
    #print line
    t=line.split()
    name2.append(t[0])
    centerv2.append(float(t[1]))
    sigma2.append(float(t[2]))
infile.close()
name2=np.array(name2,'S10')
centerv2=np.array(centerv2,'f')
sigma2=np.array(sigma2,'f')


    

In [25]:

    

centerv=np.zeros(len(name),'f')
sigma=np.zeros(len(name),'f')
for i in range(len(name)):
    centerv[i]=centerv2[name[i] == name2]
    sigma[i] = sigma2[name[i] == name2]


    

In [26]:

    

# calculate R200 for each cluster
cr200=getr200(centerv,sigma)


    

In [27]:

    

name


    

    
Out[27]:





array(['NRGb004', 'NRGs027', 'NRGs038', 'NRGs076', 'NRGs090', 'NRGs110',
       'NRGs117', 'NRGb128', 'NRGb155', 'NRGb177', 'NRGb226', 'NRGb244',
       'NRGb247', 'NRGs317', 'Abell2063'], 
      dtype='|S10')

In [32]:

    

# loop over clusters and make phase-space diagram for each
clusterswedontlike=['NRGb004','NRGs038','NRGs076','NRGs090','NRGb128']
topfive=['NRGb155','NRGb226','Abell2063','NRGs117','NRGb247']
for i in range(len(name)):
    if name[i] in clusterswedontlike:
        print "found one you don't like: ",name[i]
        continue
    if name[i] in topfive:
        cl=cluster(name[i],centerv[i],sigma[i],ra[i],dec[i],cr200[i])
        cl.readNSA()
        cl.readAGC()
        cl.plotphasespace()


    

    




found one you don't like:  NRGb004
found one you don't like:  NRGs038
found one you don't like:  NRGs076
found one you don't like:  NRGs090
found one you don't like:  NRGb128






    













    













    













    













    

In [70]:

    

if 'NRGb004' in name: print 'found it!'


    

    




found it!

In [ ]: