def drawDensityProfile(self, catalog=None):
rmax = 24. # arcmin
bins = numpy.arange(0, rmax + 1.e-10, 2.)
centers = 0.5 * (bins[1:] + bins[0:-1])
area = numpy.pi * (bins[1:]**2 - bins[0:-1]**2)
r_peak = self.kernel.extension
stars = self.get_stars()
angsep = ugali.utils.projector.angsep(self.ra, self.dec,
stars.ra, stars.dec)
angsep_arcmin = angsep * 60 # arcmin
cut_iso = self.isochrone_selection(stars)
h = numpy.histogram(angsep_arcmin[(angsep_arcmin < rmax) & cut_iso], bins=bins)[0]
h_out = numpy.histogram(angsep_arcmin[(angsep_arcmin < rmax) & (~cut_iso)], bins=bins)[0]
gals = self.get_galaxies()
if len(gals):
angsep_gal = ugali.utils.projector.angsep(self.ra, self.dec,
gals.ra, gals.dec)
angsep_gal_arcmin = angsep_gal * 60 # arcmin
cut_iso_gal = self.isochrone_selection(gals)
h_gal = np.histogram(angsep_gal_arcmin[(angsep_gal_arcmin < rmax) & cut_iso_gal], bins=bins)[0]
h_gal_out = np.histogram(angsep_gal_arcmin[(angsep_gal_arcmin < rmax) & (~cut_iso_gal)], bins=bins)[0]
plt.plot(centers, h/area, c='red', label='Filtered Stars')
plt.errorbar(centers, h/area, yerr=(numpy.sqrt(h) / area), ecolor='red', c='red')
plt.scatter(centers, h/area, edgecolor='none', c='red', zorder=22)
plt.plot(centers, h_out/area, c='gray', label='Unfiltered Stars')
plt.errorbar(centers, h_out/area, yerr=(numpy.sqrt(h_out) / area), ecolor='gray', c='gray')
plt.scatter(centers, h_out/area, edgecolor='none', c='gray', zorder=21)
if len(gals):
plt.plot(centers, h_gal/area, c='black', label='Filtered Galaxies')
plt.errorbar(centers, h_gal/area, yerr=(numpy.sqrt(h_gal) / area), ecolor='black', c='black')
plt.scatter(centers, h_gal/area, edgecolor='none', c='black', zorder=20)
plt.xlabel('Angular Separation (arcmin)')
plt.ylabel(r'Density (arcmin$^{-2}$)')
plt.xlim(0., rmax)
ymax = pylab.ylim()[1]
#pylab.ylim(0, ymax)
pylab.ylim(0, 12)
pylab.legend(loc='upper right', frameon=False, fontsize=10)
评论列表
文章目录
