Plotting routines for a set of models

Collections of plotting routines to produce comparison plots for a set of ProDiMo models.

# the prodimopy modules for reading and plotting
import prodimopy.read as pread
import prodimopy.plot_models as ppms
# this is for the PDF output
from matplotlib.backends.backend_pdf import PdfPages

# Load the prodimopy style
import matplotlib.pyplot as plt
plt.style.use("prodimopy")

# the ProDiMo model directories, each directory is one model
# adapt to your needs
# for some plotting routines it is assumed that the last model
# in that list is a reference model
mdirs=["model1","model2","model3"]

# make a list of models and read all of them
models=list()
for mdir in mdirs:
  models.append(pread.read_prodimo(mdir))

# one can also manipulate the name of a model, this name
# will be shown in the legend of the plots
models[-1].name="reference"

# Create out_models.pdf where all the various plots are stored
with PdfPages("out_models.pdf") as pdf:
  # create the PlotModels object, here it is also possible to e.g.
  # configure the colors used for line plots. For Details see the
  # documentation
  pms=ppms.PlotModels(pdf)

  pms.plot_NH(models,xlog=True)

  # if there was at least one SED calculated for the models, plot them including
  # observational data for comparison (if available)
  if not all(x.sed is None for x in models):
    pms.plot_sed(models,ylim=[1.e-17,None],plot_starSpec=False,sedObs=models[0].sedObs)


  # compare some column densities and average abundances
  for spec in ["CO","H2O"]:
    pms.plot_tcdspec(models, spec,xlog=True)
    pms.plot_avgabun(models,spec,xlog=True)


  # get all line estimates for CO in the wl range from 200 to 1500 mic
  lines=models[0].selectLineEstimates("CO",wlrange=[200,1500])


  # build an array of line idents which are considered for the plot
  # here also only CO lines with an upper level smaller 20 are considered (i.e. the main transition)
  # here jup is the internal level numbering of ProDiMo
  idents=[[x.ident,x.wl] for x in lines if x.jup <20]
  pms.plot_lines(models,idents,showBoxes=False)

class prodimopy.plot_models.PlotModels(pdf=None, colors=None, styles=None, markers=None, fs_legend=None, ncol_legend=0)[source]

Bases: object

Parameters:

pdf (PdfPages | None)
colors (list[Any] | None)
styles (list[str] | None)
markers (list[str] | None)
fs_legend (int | None)
ncol_legend (int)

pdf: A matplotlib.backends.backend_pdf.PdfPages object for output as pdf, or None, for output on screen (notebook etc.). Default: None

colors: Default colors for the different models. If not set the matplotlib default color cycle is used.

styles: Default line styles for the different models. If not set all models have solid lines (-).

markers: Default marker styles for the different models. If not set all models have D (diamond) as marker.

fs_legend: Default font size for the legend in a plot. If not set it is taken from the matplotlib rcParams.

ncol_legend: Number of entries per columns in the legend. Default: 0 (automatic)

plot_lines(models, lineIdents, useLineEstimate=True, jansky=False, showBoxes=True, lineObs=None, lineObsLabel='Obs.', peakFlux=False, incidx=0, showCont=False, xLabelGHz=False, showGrid=True, **kwargs)[source]

Plots the fluxes for a given selection of lines or lineEstimates.

Parameters:

models (list[Data_ProDiMo]) – A list of ProDiMo models.
lineIdents (list[tuple[str, float]]) – a list of line identifiers. Each entry should contain ("ident",wl) For example: lineIdents=[("CO",1300),("CO",800)]. Those values are passed to getLineEstimate() or getLine() (see next useLineEstimate). The order of the lineIdents also defines the plotting order of the lines (from left to right).
useLineEstimate (bool) – If True the lines are selected from the lineEstimates otherwise the lines are taken from the full line radiative transfer calculations. Default: True
jansky (bool) – plot the fluxes in units of Jansky. Default: True.
showBoxes (bool) – Show boxes that indicate the “deviation” region within a factor of 3 and 10 with respect to reference line flux (that list model of the list). Default: True
lineObs (list[DataLineObs] | None) – a list of DataLineObs should include all the lines that were include in the line transfer and in the same order. The observations are plotted together with the model line fluxes. Default: None (no observations are plotted)
peakFlux (bool) – If True use the peakFlux from the line profile instead of the total line flux. It can be useful for making predictions for planned observations, as in that case the peakFlux is often used to estimate observing times. If peakFlux = “convolved”, use the peak from the convolved line profile instead. Only works with useLineEstimate=False and jansky=True.
incidx (int) – Fluxes from the line ray-tracing (useLineEstimate=False) can exist for multiple inclinations. This parameter allows to chose one. The parameter is just passed on to getLine()). Default: 0
showCont (bool) – Also show the continuum fluxes at the line positions.
xLabelGHz (bool) – Use GHz in the x axis labels instead of wavelength.
showGrid (bool) – Show a grid in the plot.
lineObsLabel (str)

Todo

maybe two routines, one for lines and one for lineEstimates (sharing the commong parts)
if not FlineEstimates than it would be possible to plot all lines for which line transfer is done (would make things a bit easier)

plot_NH(models, sdscale=False, marker=None, **kwargs)[source]

Plots the total vertical hydrogen column density as a function of radius for all given models.

Parameters:

models (list(Data_ProDiMo)) – the models to plot
sdscale (boolean) – set to True to show the scale in g/cm^2 for the second y axis Default : False
marker (str) – if not None plot also the given marker symbol. Default : None

plot_tcdspec(models, species, relToH=False, facGrayBox=3.0, ice=False, linewidths=None, ax=None, **kwargs)[source]

Plots the total vertical column density for the given species for all the models as a function of the radius.

Parameters:

models (list(Data_ProDiMo)) – the models to plot
species (str) – the name of the species to plot.
relToH (boolean) – Plot the columnd density relative to the total hydrogen column density. Which is essential equal to the average abundance at a certain radius. DEFAULT: False
facGrayBox (float) – Plot a gray area for each column density line with width given by y / facGrayBox and y * facGrayBox. Uses matplotlib:fill_between. Set to None if it should not be plotted. DEFAULT: 3
ice (boolean) – Plot also the ice phase column density if the ice species exists.
linewidths (array_like) – An array of linewidths for each model. If None the default linewidth is used. If linewidths is set it has to have the same lengths as the models array.
ax (Axes) – If not None the plot is made on the given Axes object, otherwise a new Figure is created.

plot_tauline(models, lineIdent, xlog=True, **kwargs)[source]

Plots the line optical depths as a function of radius for a given line.

Parameters:

models (list(Data_ProDiMo)) – the models to plot
lineIdent (array_like(ndim=1)) – two element array. First item lineIdent (e.g. CO) second item wavelength (mu m) see also getLineEstimate()
FIXME (better treatment if no lines are found (should not crash).)

plot_avgabun(models, species, **kwargs)[source]

Plots the average abundance of the species as a function of radius. The average abundance is given by the vertical column density of the species divided by the total hydrogen column density.

Parameters:

models (list(Data_ProDiMo)) – The models to plot.
species (str) – The name of the species to plot.

plot_abunvert(models, r, species, **kwargs)[source]

Plot abundances at the given radius r as a function of height.

Parameters:

models (list(Data_ProDiMo)) – The models to plot.
r (float) – The radius at which the plot should be made UNIT: [au]
species (str) – The species name for which the abundance should be plotted.

Todo

Make this a wrapper of plot_vertical(). Just pass nmol as fieldname and species name

plot_radial(models, fields, ylabel=None, zidx=0, showlegend=True, ax=None, **kwargs)[source]

Plots a quantity in radial direction. Fields must have the same number of entries as models and must contain arrays with the dimension of nx.

Parameters:

models (list[Data_ProDiMo]) – the models to plot.
fields (list[ndarray] | str) –
If a string is passed, it is interpreted as the name of attribute (e.g. "td", "rhog") of Data_ProDiMo.

If passed as list, each entry in the list needs to be an array with dimension of NX and corresponds to one model. The order has to be the same as for the models list.

Example:.
```
zidx=5
# Using a attribute names
ppm.plot_radial(models,"rhog",zidx=zidx)
# Using a fields array, this allows also to manipulate the values, if that is usefull.
fields=[mod.rhog[:,zidx] for mod in models]
ppm.plot_radial(models,fields,ylabel="rho [g/cm^3]")
```
ylabel (str | None)
zidx (int)
showlegend (bool)
ax (Axes | None)

Return type:

Figure | None

ylabel :: The label for the y Axis of the plot. if fields is a string, any ylabel is None, an automatic label is generated.
zidx :: The z index for which height the radial plot should be made. This is used to properly calculate the radius. DEFAULT: z=0
showlegend :: If True a legend is shown. DEFAULT: True

ax : Use the passed axes to make the plot. If None a new Figure is created.

Todo

implement things from plot_midplane (e.g. species)

plot_midplane(models, fieldname=None, ylabel=None, xfieldname=None, xlabel=None, species=None, patches=None, scaleToRin=False, ax=None, **kwargs)[source]

Plots a quantity in in the midplane as a function of radius.

Parameters:

models (list(Data_ProDiMo)) – the models to plot
fieldname (str | None) – fieldname is any attribute of Data_ProDiMo. If None the species parameter has to be set.
ylabel (str | None) – The label for the y axis of the plot.
xfieldname (str | None) – Is the name of any attribute of Data_ProDiMo and will be used as the x coordinate in the plot. If xfieldName=None the x coordinate of the models is used. DEFAULT: None
xlabel (str | None) – The label for the x-coordinate of the plot.
species (str | None) – The name of a chemical species. If set, the abundance of this species is plotted. DEFAULT: None
patches (matplotlib.patches object) – A matplotlib patches object that is overplotted on the figure.
scaleToRin (bool) – If True the r-Rin is used for the x-coordinate of the plot.
ax (Axes | None) – Used the passed axes to make the plot. If None a new Figure is created.
todo:: (..) –
- merge with plot_radial()
- or alternatively call plot_radial with the proper parameters

Return type:

Figure | None

plot_vertical_nH(models, r, field, ylabel, species=None, patches=None, showR=True, **kwargs)[source]: Deprecated since version 2.5.1: Use plot_vertical() with xfield=’NHver’ instead.

plot_vertical(models, r, fieldname=None, ylabel=None, fieldnameIdx=None, species=None, xfield='zr', ylog=True, xlog=False, zr=True, showlegend=True, ax=None, **kwargs)[source]

Plots a quantity (fieldname) as a function of height (various options for the coordinates) at the given radius.

Parameters:

models (list[Data_ProDiMo]) – The list of models to plot.
r (float) – The radius at which the vertical cut should be made (uses nearest neighbour). Unit: au.
fieldname (str | None) – The name of a field(attribute) of the model data object (Data_ProDiMo). If not set species has to be provided.
fieldnameIdx (int | None) – An additional index in case the quantity to plot has more than the two spatial dimensions. e.g. field[ix,iz,fieldnameIdx]
species (str | None) – A species name. If provided the abundance of this species is plotted. fieldname has to be nmol in that case.
xfield (str) –
Chose the field/quantity that should be used for the x-axis (height). Options:
- zr z/r, default
- NHver vertical hydrogen column density.
- AVver vertical visual extinction.
- tg gas temperature.
If none of those or an unknown value is given the z coordinate is used.
ylabel (str | None)
ylog (bool)
xlog (bool)
zr (bool)
showlegend (bool)
ax (Axes | None)

Return type:

Figure | None

plot_abunradial(models, species, perH2=False, **kwargs)[source]

Plots the abundance of the given species as a function of radius for the midplane of the disk.

Parameters:

models (list(Data_ProDiMo)) – the models to plot
species (str) – The name of the species to plot.
perH2 (boolean) – Plot the abundance relative to molecular hydrogen number density instead of the total hydrogen number density.

plot_sed(models, plot_starSpec=True, sedObs=None, sedObsModels=False, unit='erg', incidx=0, reddening=False, ax=None, **kwargs)[source]

Plots the Seds and the StarSpectrum (optionally).

Parameters:

models (list(Data_ProDiMo)) – the models to plot
plot_starSpec (boolean) – If True also plot the stellar spectrum.
sedObs (DataContinuumObs) – One set of the SED observations to plot (photometry)
sedObsModels (boolean) – If True use the DataContinuumObs object from each model individually. So for each model the provided observations are plotted.
incidx (int) – which inclination (index) should be used for plotting (0 is the first one and default). Currently only a single value is supported, as otherwise the plotting becomes messy.
unit (str) – In which unit the SED should be plotted. Allowed values: erg [erg/s/cm^2], W [W/m^2], Jy, mJy.
reddening (boolean or str) – If True also apply the reddening of the SED using the given A_V value from the provided observational data. If reddening is a string the according extinction model is used in extinction().
ax (matplotlib.axes) – The axes object to plot the SEDs on. Default: None (new figure is created)

plot_dust_opac(models, xenergy=False, opactype='both', pseudoaniso=False, **kwargs)[source]

Plots the dust opacities, absorption, scattering or extinction.

Parameters:

opactype (str) – both absorption and scattering (separately), abs only absorption, sca only scattering, ext only extinction
pseudoaniso (boolean) – Use the pseudo anisotropic scattering opacity (Default: False).

plot_starspec_xray(models, nuInu=False, ax=None, **kwargs)[source]: Plots the X-ray spectrum.

plot_starspec(models, step=10, xunit='micron', nuInu=True, linewidth=1.0, ax=None, second_xaxis=False, **kwargs)[source]

Plots the full Stellar Spectrum for each model.

Parameters:

models (list[Data_ProDiMo]) – the models to plot
step (int) – Only plot every step point of the spectrum.
xunit (str) – Unit for the x-axis. Currently mircon and eV are possible.
nuInu (bool) – Plot \(\nu\) times \(I_\nu\) instead of \(I_\nu\) only.
linewidth (float) – Line width for the spectrum lines.
ax (Axes) – The axes object to plot the spectra on. Default: None (new figure)
second_xaxis (bool) – If True add a second x-axis showing the corresponding energy/wl in eV/micron.

plot_line_profil(models, wl, ident=None, linetxt=None, lineObs=None, unit='Jy', normalized=False, **kwargs)[source]: Deprecated since version 2.4.0: Use plot_lineprofile() instead.

plot_lineprofile(models, wl, ident=None, contsub=True, linetxt=None, lineObs=None, incidx=0, unit=None, normalized=False, convolved=False, style=None, ax=None, **kwargs)[source]

Plots the line profile for the given line (id wavelength and line ident)

Parameters:

models (array_like(Data_ProDiMo,dim=1)) – the models
wl (float) – The wavelength of the line in micrometer. Plotted is the line closest to wl.
ident (str) – The optional line ident which is additionally use to identify the line.
contsub (boolean) – Remove the continuum. Default: True
linetxt (str) – A string the is used as the label for the line. Default: None
lineObs (array_like(ndim=1)) – list of prodimopy.read.DataLineObs objects. Must be consistent with the list of lines from the line radiative transfer.
incidx (int) –
which inclination (index) should be used for plotting (0 is the first one and default).

If lineObj is passed only one inclination (the one set in that obj) will be plotted.
unit (str | None) – In what unit should the line flux be plotted (see flux_unit). if None the current unit of the line profile is used.
normalized (boolean) – if True normalize the profile to the peak flux of each line
convolved (boolean) – if True plot the convolved profile.
style (str) – if style is step the profile is plotted as a step function assuming the values are the mid point of the bin.

prodimopy.plot_models.load_style(style='prodimopy')[source]

Simple wrapper that calls load_mplstyle() with the given style.

Parameters:: style (str) – The name of the style to load. Default: “prodimopy”