Dressler Shectman

Information taken from: "Evidence for Substructure in Rich Clusters of Galaxies from Radial-Velocity Measurements" by Alan Dressler and Stephen Shectman

For each galaxy with a measured radial velocity, we find the ten nearest neighbors with measured velocities. From this sample of 11, we calculate the local velocity mean and dispersion determined for the entire cluster sample. The global mean velocity $\vec{v}$ and dispersion $\sigma$ are calculated assuming a Gaussian distribution, but this assumption is ultimately removed by the normalizing procedure.

$$\delta^2 = (11/\sigma^2)[(\bar{v_{local}} - \bar{v})^2 + (\sigma_{local} - \sigma)^2]$$

In [1]:
import glob
from matplotlib import pyplot as plt
import numpy as np
import scipy.stats
import glob
from astropy.stats import biweight_location, biweight_midvariance, sigma_clip, bootstrap
from astropy import cosmology
from astropy.io import fits
from math import sqrt

%matplotlib inline


/Applications/anaconda/lib/python2.7/site-packages/IPython/kernel/__init__.py:13: ShimWarning: The `IPython.kernel` package has been deprecated. You should import from ipykernel or jupyter_client instead.
  "You should import from ipykernel or jupyter_client instead.", ShimWarning)

In [28]:
def distance(ra,dec,all_ra,all_dec, n=10):
    d = np.sqrt((ra-all_ra)**2+(dec-all_dec)**2)
    t = np.argsort(d)
    return t[0:n+1]

if __name__ == '__main__':

    cl = fits.getdata('Abell2063_NSA.fits')

    # cols = A2063_AGC[1].columns # Prints the titles of each column
    # cols.info()

    ds = np.zeros(len(cl.RA))

    
    ave_v=np.mean(cl.ZDIST*3.e5)
    std_v=np.std(cl.ZDIST*3.e5)
    plt.figure()
    plt.plot(cl.RA,cl.DEC,'k.')
    #for i in range(1):
    

    for i in range(len(cl.RA)):
        d_index = distance(cl.RA[i],cl.DEC[i],cl.RA,cl.DEC, n=10)
        #print d_index

        # calculate avg velocity
        ave_local = np.mean(cl.ZDIST[d_index]*3.e5)
        std_local = np.std(cl.ZDIST[d_index]*3.e5)
        
        ds[i] = 11./std_v**2*((ave_local - ave_v)**2+ (std_local - std_v)**2)


        if i == 250:
            plt.plot(cl.RA[d_index],cl.DEC[d_index],'ro',markersize=10)
            plt.plot(cl.RA[i],cl.DEC[i],'bs')
            print ds[i]
            plt.scatter(cl.RA[i],cl.DEC[i],s=ds[i],color='k') 

    
plt.figure()
plt.scatter(cl.RA,cl.DEC,s=ds,c=cl.ZDIST*3.e5)
plt.colorbar()
#plt.axis([230,232,7,10])


0.0396990532493
Out[28]:
<matplotlib.colorbar.Colorbar instance at 0x10fdc2950>