Use center values for interpolation

This commit is contained in:
Noe Brucy
2019-05-14 16:57:16 +02:00
parent 18cfbc6eb3
commit a587704ed5
+17 -20
View File
@@ -675,7 +675,7 @@ def disk_prop(
"time": time, "time": time,
"mass_disk": total_mass_disk, "mass_disk": total_mass_disk,
"mass_box": total_mass, "mass_box": total_mass,
"rad": rad[0 : nb_bin - 1], "rad": rad[:-1],
"center": pos_star, "center": pos_star,
"alpha_rey": alpha_rey_rad, "alpha_rey": alpha_rey_rad,
"alpha_rey_mean": alpha_rey_mean, "alpha_rey_mean": alpha_rey_mean,
@@ -873,8 +873,7 @@ def disk_pdf(
nb_bin_hist=50, nb_bin_hist=50,
tag="", tag="",
rad_min=0.08, rad_min=0.08,
rad_max=0.35, rad_max=0.4,
interpol=True,
): ):
# Load property file # Load property file
@@ -905,8 +904,18 @@ def disk_pdf(
coldens_map = maps_disk["coldens_z"] coldens_map = maps_disk["coldens_z"]
im_extent = maps_disk["im_extent"] im_extent = maps_disk["im_extent"]
# radius of the corner of the box plus a margin
rad_box = (
np.sqrt((im_extent[1] - pos_star[0]) ** 2 + (im_extent[3] - pos_star[1]) ** 2)
+ 0.1
)
# radial bins
rad_bins = prop_disk["rad"] rad_bins = prop_disk["rad"]
rad_bins = 0.5 * (rad_bins[0:-1] + rad_bins[1:])
rad_bins = np.concatenate(([0.0], rad_bins, [rad_box]))
coldens_rad = prop_disk["coldens"] coldens_rad = prop_disk["coldens"]
# Add value for borders
coldens_rad = np.concatenate(([coldens_rad[0]], coldens_rad))
x = np.linspace(im_extent[0], im_extent[1], coldens_map.shape[0]) x = np.linspace(im_extent[0], im_extent[1], coldens_map.shape[0])
y = np.linspace(im_extent[2], im_extent[3], coldens_map.shape[1]) y = np.linspace(im_extent[2], im_extent[3], coldens_map.shape[1])
@@ -920,25 +929,13 @@ def disk_pdf(
bins_flat = bins.flatten() bins_flat = bins.flatten()
rr_flat = rr.flatten() rr_flat = rr.flatten()
# Retrieve the mean of the bi # Retrieve the mean of the bins
if interpol: # use linear interpolation to improve accuracy # use linear interpolation to improve accuracy
rad_bins_ext = np.append( coldens_mean = (rad_bins[bins_flat + 1] - rr_flat) * coldens_rad[bins_flat]
rad_bins,
np.sqrt(
(im_extent[1] - pos_star[0]) ** 2 + (im_extent[3] - pos_star[1]) ** 2
),
)
coldens_rad_ext = np.append(coldens_rad, coldens_rad[-1])
coldens_mean = (rad_bins_ext[bins_flat + 1] - rr_flat) * coldens_rad[bins_flat]
coldens_mean = ( coldens_mean = (
coldens_mean coldens_mean + (rr_flat - rad_bins[bins_flat]) * coldens_rad[bins_flat + 1]
+ (rr_flat - rad_bins[bins_flat]) * coldens_rad_ext[bins_flat + 1]
) )
coldens_mean = coldens_mean / ( coldens_mean = coldens_mean / (rad_bins[bins_flat + 1] - rad_bins[bins_flat])
rad_bins_ext[bins_flat + 1] - rad_bins[bins_flat]
)
else:
coldens_mean = coldens_rad[bins_flat]
coldens_mean_map = np.reshape(coldens_mean, coldens_map.shape) coldens_mean_map = np.reshape(coldens_mean, coldens_map.shape)