Plot refactor#415
Conversation
| cmap = colors.LinearSegmentedColormap.from_list( | ||
| "trunc({n},{a:.2f},{b:.2f})".format(n=cmap.name, a=0.0, b=0.7), | ||
| cmap(np.linspace(0.0, 0.7, 1000)), | ||
| ) |
There was a problem hiding this comment.
Do we need to define cmap here and below? Can we reduce it to just one def, either at the top here or in the Plotter class?
| def __init__( | ||
| self, | ||
| results, | ||
| object_to_plot=1, # TODO: Support multiplanet |
There was a problem hiding this comment.
Noting that we should resolve this todo before merge (I'm sure this is already on your radar).
| self.results = results | ||
| self.system = results.system | ||
|
|
||
| if default_cmap is not None: |
There was a problem hiding this comment.
Make sure the logic works here to define the colormap as we expect
| self.cmap = default_cmap | ||
|
|
||
| if start is None: | ||
| start = getattr(Time(np.min(self.system.data_table['epoch'])-365*3, format="mjd"), time_format) |
There was a problem hiding this comment.
I'm tempted to automatically set the start and end times for each dataset being plotted, and asking the user to manually override these after the Figure object is returned. If we make this change, we should include an example of manually setting axis limits in the plotting tutorial. Also, if there is no data, but the user still wants to plot the prediction, we can use the "first epoch (of any kind) - 3 years" as the start date, and analagous for end date.
| cbar_param=None, | ||
| rv_time_series=None, | ||
| rv_time_series2=None, | ||
| ): |
| ) | ||
|
|
||
|
|
||
| def legacy_plot_orbits( |
There was a problem hiding this comment.
Decided you can remove the legacy code. It will be in the git history if we want to refer back to it.
| fig=None, | ||
| ): | ||
| # TODO: Improve message | ||
| print("The orbitize.plot.plot_orbits function is deprecated. Instead, initialize an orbitize.plot.Plotter and use the plot_orbits method.") |
There was a problem hiding this comment.
Very nice. Just writing down what we talked about over zoom for posterity that we should raise a DeprecationWarning (or similar) here
| self.sep_data, self.sep_err, self.pa_data, self.pa_err, self.ra_data, self.ra_err, self.dec_data, self.dec_err = self._calc_seppa_radec(self.data) | ||
| self.astr_raoff, self.astr_deoff, self.astr_vz, self.astr_epochs = self._calc_astr_orbits(self.standard_post, self.num_orbits_to_plot, self.object_to_plot) | ||
|
|
||
| def _calc_seppa_radec(self, all_data): |
There was a problem hiding this comment.
add docstrings for all these private functions
| self.astr_raoff, self.astr_deoff, self.astr_vz, self.astr_epochs = self._calc_astr_orbits(self.standard_post, self.num_orbits_to_plot, self.object_to_plot) | ||
|
|
||
| def _calc_seppa_radec(self, all_data): | ||
| data = all_data[all_data["quant_type"]!="rv"] |
There was a problem hiding this comment.
can you calculate all of these data type masks outside of these private methods? It would be nice to see them in one place (maybe in set_params above)
There was a problem hiding this comment.
Actually, I disagree with myself now that I went through the whole refactor. I like it as is.
| # TODO: vectorize | ||
| for i in np.arange(num_orbits_to_plot): | ||
| # Calculate ra/dec offsets for all epochs of this orbit | ||
| raoff0, deoff0, vz0 = self.system.compute_all_orbits( |
There was a problem hiding this comment.
note we just changed the return type of compute_all_orbits in v3.4.0, so you will need to add raoff0, deoff0, vz0, _ here (and elsewhere this method is called)
|
|
||
| return raoff, deoff, vz, astr_epochs | ||
|
|
||
| def _calc_rv(self, object): |
There was a problem hiding this comment.
we should change object to something else since object is a special variable in python. It should work fine as is but might be confusing.
| return rv_data, insts, gams, labels, gam_idx, inds, sig_idx | ||
|
|
||
| def _get_standard_post(self, num_orbits_to_plot): | ||
| # TODO: Replace random with results.downsample |
There was a problem hiding this comment.
good call to replace with results.downsample!
| plot_astrometry_insts=False, | ||
| rv_time_series=False, | ||
| rv_time_series2=False, | ||
| # plot_errorbars=True, |
| warnings.simplefilter("ignore", ErfaWarning) | ||
|
|
||
| # Before starting to plot rv data, make sure rv data exists: | ||
| rv_indices = np.where(self.data["quant_type"] == "rv") |
There was a problem hiding this comment.
Can we replace this bit of code with the self.rv_data variable you defined above?
There was a problem hiding this comment.
Actually, I think we should raise this warning in self.set_params.
| else: | ||
| ax = plt.subplot2grid((2, 14), (0, 0), rowspan=2, colspan=6) | ||
|
|
||
| astr_inds = np.where((~np.isnan(self.data["quant1"])) & (~np.isnan(self.data["quant2"]))) |
There was a problem hiding this comment.
Let's use the variables you defined above rather than recalculating this mapping here.
|
|
||
| # For plotting different astrometry instruments | ||
| if plot_astrometry_insts: | ||
| astr_colors = ("#FF7F11", "#11FFE3", "#14FF11", "#7A11FF", "#FF1919") |
There was a problem hiding this comment.
define these in the same place you define the RV data colors and symbols above.
|
|
||
| return fig | ||
|
|
||
| def _plot_full_orbits(self, ax, plot_astrometry, square_plot, fontsize, cmap, plot_astrometry_insts, astr_insts=None, astr_inst_inds=None, astr_colors=None, astr_symbols=None): |
There was a problem hiding this comment.
I'd suggest putting all the private methods together in the beginning
| plt.sca(ax3) | ||
|
|
||
| # scale back to primary RV semi amplitude | ||
| vz0 = self.vz1[i] * (-(mtot[i] - m0[i]) / np.median(m0[i])) # TODO: vectorize |
There was a problem hiding this comment.
we should use the output of calc_all_orbits here rather than doing the math ourselves
| else: | ||
| plt.sca(ax3) | ||
|
|
||
| # scale back to primary RV semi amplitude |
There was a problem hiding this comment.
Looks like this comment isn't relevant anymore; we can remove it
| plt.plot(yr_epochs, pas, color=sep_pa_color) | ||
|
|
||
| def _plot_rv_data(self, ax3, ax4, rv_time_series, rv_time_series2, sep_pa_color): | ||
| # plot RV orbits here |
There was a problem hiding this comment.
let's call these _plot_rv_model, _plot_sep_pa_model, etc, since the data is being plotted in a separate place
| epochs = inst_data["epoch"] | ||
| epochs = Time(epochs, format="mjd").decimalyear | ||
| inst_err = inst_data["quant1_err"] | ||
| gam_err2 = np.sqrt(np.square(inst_err)+np.square(med_sigma[i])) |
There was a problem hiding this comment.
these (gam_err1, gam_err2, gam_err_2) could use more informative variable names
| s=30, | ||
| marker=next(ax3_symbols), | ||
| c=c, | ||
| label=name, |
There was a problem hiding this comment.
if the name is "defrv" (our default when the user doesn't provide an instrument name) we should make the label something more informative like "RV data"
| epochs = Time(epochs, format="mjd").decimalyear | ||
| inst_err = inst_data["quant1_err"] | ||
| gam_err2 = np.sqrt(np.square(inst_err)+np.square(med_sigma[i])) | ||
| gam_err_1 = np.sqrt(gam_err2) |
There was a problem hiding this comment.
subtle point here-- each orbit may have its own fitted value of jitter (sigma). We want to use that value, rather than the median across all orbits.
There was a problem hiding this comment.
Also, I think gam_err_1 is unneeded. Below, gam_err_2 should be np.sqrt(gam_err2**2 + stddev_ga[i]**2)
| plt.errorbar( | ||
| x=epochs2, | ||
| y=rvs2-med_ga2[i], | ||
| yerr=stddev_ga2[i], |
There was a problem hiding this comment.
Nice job working through all the differences here!
| ecolor=c, | ||
| zorder=5, | ||
| ls="none", | ||
| label="{} gamma err".format(name2) |
There was a problem hiding this comment.
I'd probably suggest writing these captions out in latex to make them as clear as possible: std(rv_{inst name} - gamma_{inst name}); sigma_{inst name}2 + sigma2 (and same for the other case below)
| deoff = self.astr_deoff | ||
| seps = [] | ||
| pas = [] | ||
| raoff_100 = self.raoff1 |
There was a problem hiding this comment.
use your defined variables throughout rather than redefining them (for clarity)
|
|
||
| # for j in range(len(astr_epochs)): | ||
|
|
||
| # seps0, pas0 = orbitize.system.radec2seppa(raoff[i][j], deoff[i][j], mod180=mod180) |
| arcs with PAs that cross 360 deg during observations (default: False) | ||
|
|
||
| Return: | ||
| ``matplotlib.pyplot.Figure``: the residual plots |
There was a problem hiding this comment.
Let's use the same logic for plotting residuals as we used above for plotting the gammas for RVs. Offer two options for plotting the uncertainties, exactly the same way you did above.
| timestep = epoch_in_yr[1] - epoch_in_yr[0] | ||
| # dra/dt and ddec/dt | ||
| ddec_b = np.gradient(deoff[i, :], timestep) # in mas/yr | ||
| dec_b_radian = ( |
There was a problem hiding this comment.
let's use astropy.units to do these conversions for clarity
|
Looks great! Excellent job @eshelDror. A few small overall comments:
|
|
|
||
| if plot_priors: | ||
| axes = figure.axes | ||
| num_params = len(param_indices) |
There was a problem hiding this comment.
I think you can avoid much of this normalization code by passing hist2d_kwargs = {"density":True} into corner.corner.
| for i, param_i in enumerate(param_indices): | ||
| prior = results.system.sys_priors[param_i] | ||
| if not hasattr(prior, "compute_lnprob"): | ||
| continue |
There was a problem hiding this comment.
This should never happen; Prior objects by definition (by means of the abstract base class) must have a compute_lnprob method, so I think you can safely remove this check
|
Also, see aesthetic suggestions in #378 and please implement the ones you like! I think the main things I'd like to see here based on those suggestions are:
|
No description provided.