def observableFractionCDF(self, mask, distance_modulus, mass_min=0.1):
"""
Compute observable fraction of stars with masses greater than mass_min in each
pixel in the interior region of the mask. Incorporates simplistic
photometric errors.
ADW: Careful, this function is fragile! The selection here should
be the same as mask.restrictCatalogToObservable space. However,
for technical reasons it is faster to do the calculation with
broadcasting here.
ADW: This function is currently a rate-limiting step in the likelihood
calculation. Could it be faster?
"""
method = 'step'
mass_init,mass_pdf,mass_act,mag_1,mag_2 = self.sample(mass_min=mass_min,full_data_range=False)
mag_1 = mag_1+distance_modulus
mag_2 = mag_2+distance_modulus
mask_1,mask_2 = mask.mask_roi_unique.T
mag_err_1 = mask.photo_err_1(mask_1[:,np.newaxis]-mag_1)
mag_err_2 = mask.photo_err_2(mask_2[:,np.newaxis]-mag_2)
# "upper" bound set by maglim
delta_hi_1 = (mask_1[:,np.newaxis]-mag_1)/mag_err_1
delta_hi_2 = (mask_2[:,np.newaxis]-mag_2)/mag_err_2
# "lower" bound set by bins_mag (maglim shouldn't be 0)
delta_lo_1 = (mask.roi.bins_mag[0]-mag_1)/mag_err_1
delta_lo_2 = (mask.roi.bins_mag[0]-mag_2)/mag_err_2
cdf_1 = norm_cdf(delta_hi_1) - norm_cdf(delta_lo_1)
cdf_2 = norm_cdf(delta_hi_2) - norm_cdf(delta_lo_2)
cdf = cdf_1*cdf_2
if method is None or method == 'none':
comp_cdf = cdf
elif self.band_1_detection == True:
comp = mask.mask_1.completeness(mag_1, method=method)
comp_cdf = comp*cdf
elif self.band_1_detection == False:
comp =mask.mask_2.completeness(mag_2, method=method)
comp_cdf = comp*cdf
else:
comp_1 = mask.mask_1.completeness(mag_1, method=method)
comp_2 = mask.mask_2.completeness(mag_2, method=method)
comp_cdf = comp_1*comp_2*cdf
observable_fraction = (mass_pdf[np.newaxis]*comp_cdf).sum(axis=-1)
return observable_fraction[mask.mask_roi_digi[mask.roi.pixel_interior_cut]]
评论列表
文章目录
