The goal of this program is to match the AGC galaxies to the NSA catalogs.
Ideally, there will be a one-to-one match, but of course, nothing is ideal...
So what do we do with galaxies that don't have an AGC match?
not sure yet...
In [1]:
from matplotlib import pyplot as plt
import numpy as np
from numpy.lib.recfunctions import append_fields
from astropy.io import fits
import glob
%matplotlib inline
In [8]:
def findnearest(x1,y1,x2,y2,delta):#use where command
matchflag=1
nmatch=0
d=np.sqrt((x1-x2)**2 + (y1-y2)**2)#x2 and y2 are arrays
index=np.arange(len(d))
t=index[d<delta]
matches=t
if len(matches) > 0:
nmatch=len(matches)
if nmatch > 1:
imatch=index[(d == min(d[t]))][0]
else:
imatch=matches[0]
else:
imatch = 0
matchflag = 0
return imatch, matchflag,nmatch
In [3]:
def compare_samples(x1,y1,x2,y2):
plt.figure()
plt.plot(x1,y1,'bo')
plt.plot(x2,y2,'ks')
In [10]:
# match radiuss = 3 arcsec
matchradius=3./3600 # convert to degrees
clusterfiles=glob.glob('*_NSA.fits')
for cfile in clusterfiles:
nsadat=fits.getdata(cfile)
name=cfile.split('_')[0]
agcfile=name+'_AGC.fits'
print 'running cluster ',name
print agcfile
agcdat=fits.getdata(agcfile)
imatch=np.zeros(len(nsadat.RA),'i')
matchflag=np.zeros(len(nsadat.RA),'bool')
nmatch=np.zeros(len(nsadat.RA),'i')
for i in range(len(nsadat.RA)):
imatch[i],matchflag[i],nmatch[i] = findnearest(nsadat.RA[i],nsadat.DEC[i],agcdat.RA,agcdat.DEC,matchradius)
outfile=name+'_matchAGCtoNSA.fits'
matchedarray=np.zeros(len(agcdat),dtype=agcdat.dtype)
matchedarray[matchflag] = agcdat[imatch[matchflag]]
#t=append_fields(testarray,'agcmatchflag',data=matchflag)
fits.writeto(outfile,matchedarray,clobber=True)
In [60]:
In [5]:
plt.figure()
plt.hist(agcdat.VOPT,color='b',bins=np.arange(0,12000,500),label='AGC')
plt.hist(nsadat.ZDIST*3.e5,color='g',bins=np.arange(0,12000,500),label='NSA')
plt.xlabel('Recession Velocity')
Out[5]:
In [ ]:
In [ ]: