def drawMembership(self, ax=None, radius=None, zidx=0, mc_source_id=1):
if not ax: ax = plt.gca()
import ugali.analysis.scan
filename = self.config.mergefile
logger.debug("Opening %s..."%filename)
f = pyfits.open(filename)
distance_modulus = f[2].data['DISTANCE_MODULUS'][zidx]
for ii, name in enumerate(self.config.params['isochrone']['infiles']):
logger.info('%s %s'%(ii, name))
isochrone = ugali.isochrone.Isochrone(self.config, name)
mag = isochrone.mag + distance_modulus
ax.scatter(isochrone.color,mag, color='0.5', s=800, zorder=0)
pix = ang2pix(self.nside, self.glon, self.glat)
likelihood_pix = ugali.utils.skymap.superpixel(pix,self.nside,self.config.params['coords']['nside_likelihood'])
config = self.config
scan = ugali.analysis.scan.Scan(self.config,likelihood_pix)
likelihood = scan.likelihood
distance_modulus_array = [self.config.params['scan']['distance_modulus_array'][zidx]]
likelihood.precomputeGridSearch(distance_modulus_array)
likelihood.gridSearch()
p = likelihood.membershipGridSearch()
sep = ugali.utils.projector.angsep(self.glon, self.glat, likelihood.catalog.lon, likelihood.catalog.lat)
radius = self.radius if radius is None else radius
cut = (sep < radius)
catalog = likelihood.catalog.applyCut(cut)
p = p[cut]
cut_mc_source_id = (catalog.mc_source_id == mc_source_id)
ax.scatter(catalog.color[cut_mc_source_id], catalog.mag[cut_mc_source_id], c='gray', s=100, edgecolors='none')
sc = ax.scatter(catalog.color, catalog.mag, c=p, edgecolors='none')
ax.set_xlim(likelihood.roi.bins_color[0], likelihood.roi.bins_color[-1])
ax.set_ylim(likelihood.roi.bins_mag[-1], likelihood.roi.bins_mag[0])
ax.set_xlabel('Color (mag)')
ax.set_ylabel('Magnitude (mag)')
try: ax.cax.colorbar(sc)
except: pylab.colorbar(sc)
评论列表
文章目录