Plotting

gas

Functions for plotting gas quantities

pynbody.plot.gas.rho_T(sim, rho_units=None, rho_range=None, t_range=None, **kwargs)[source]

Plot Temperature vs. Density for the gas particles in the snapshot.

Optional keywords:

rho_units: specify the density units (default is the same units as the current ‘rho’ array)

t_range: list, array, or tuple
size(t_range) must be 2. Specifies the temperature range.
rho_range: tuple
size(rho_range) must be 2. Specifies the density range.

See hist2d() for other plotting keyword options

pynbody.plot.gas.temp_profile(sim, center=True, r_units='kpc', bin_spacing='equaln', clear=True, filename=None, **kwargs)[source]

Centre on potential minimum, align so that the disk is in the x-y plane, then plot the temperature profile as a function of radius.

generic

Flexible and general plotting functions

pynbody.plot.generic.hist2d(xo, yo, weights=None, mass=None, gridsize=(100, 100), nbins=None, make_plot=True, **kwargs)[source]

Plot 2D histogram for arbitrary arrays that get passed in.

Input:

x: array

y: array

Optional keywords:

x_range: list, array, or tuple
size(x_range) must be 2. Specifies the X range.
y_range: tuple
size(y_range) must be 2. Specifies the Y range.
gridsize: (int, int) (default (100,100))
number of bins to use for the 2D histogram
nbins: int
number of bins for the histogram - if specified, gridsize is set to (nbins,nbins)
nlevels: int
number of levels to use for the contours
logscale: boolean
whether to use log or linear spaced contours
weights: numpy array of same length as x and y
if weights is passed, color corresponds to the mean value of weights in each cell
mass: numpy array of masses same length as x andy
must also have weights passed in. If you just want to weight by mass, pass the masses to weights
colorbar: boolean
draw a colorbar
scalemin: float
minimum value to use for the color scale
scalemax: float
maximum value to use for the color scale
pynbody.plot.generic.gauss_kde(xo, yo, weights=None, mass=None, gridsize=(100, 100), nbins=None, make_plot=True, nmin=None, nmax=None, **kwargs)[source]

Plot 2D gaussian kernel density estimate (KDE) given values at points (x, y).

Behavior changes depending on which keywords are passed:

If a weights array is supplied, produce a weighted KDE.

If a mass array is supplied, a mass density is computed.

If both weights and mass are supplied, a mass-averaged KDE of the weights is computed.

By default, norm=False is passed to fast_kde() meaning that the result returned is not normalized such that the integral over the area equals one.

Since this function produces a density estimate, the units of the

output grid are different than for the output of hist2d(). To get to the same units, you must multiply by the size of the cells.

Input:

xo: array

yo: array

Optional keywords:

mass: numpy array of same length as x and y
particle masses to be used for weighting
weights: numpy array of same length as x and y
if weights is passed, color corresponds to the mean value of weights in each cell
nmin: float (default None)
if weights and mass are both specified, the mass-weighted contours are only drawn where the mass exceeds nmin.
gridsize: (int, int) (default: 100,100)
size of grid for computing the density estimate
nbins: int
number of bins for the histogram - if specified, gridsize is set to (nbins,nbins)
make_plot: boolean (default: True)
whether or not to produce a plot

Keywords passed to fast_kde():

norm: boolean (default: False)
If False, the output is only corrected for the kernel. If True, the result is normalized such that the integral over the area yields 1.
nocorrelation: (default: False) If True, the correlation
between the x and y coords will be ignored when preforming the KDE.

Keywords passed to make_contour_plot():

x_range: list, array, or tuple
size(x_range) must be 2. Specifies the X range.
y_range: tuple
size(y_range) must be 2. Specifies the Y range.
nlevels: int
number of levels to use for the contours
logscale: boolean
whether to use log or linear spaced contours
colorbar: boolean
draw a colorbar
scalemin: float
minimum value to use for the color scale
scalemax: float
maximum value to use for the color scale
pynbody.plot.generic.make_contour_plot(arr, xs, ys, x_range=None, y_range=None, nlevels=20, logscale=True, xlogrange=False, ylogrange=False, subplot=False, colorbar=False, ret_im=False, cmap=None, clear=True, legend=False, scalemin=None, levels=None, scalemax=None, filename=None, **kwargs)[source]

Plot a contour plot of grid arr corresponding to bin centers specified by xs and ys. Labels the axes and colobar with proper units taken from x

Called by hist2d() and gauss_density().

Input:

arr: 2D array to plot

xs: x-coordinates of bins

ys: y-coordinates of bins

Optional Keywords:

x_range: list, array, or tuple (default = None)
size(x_range) must be 2. Specifies the X range.
y_range: tuple (default = None)
size(y_range) must be 2. Specifies the Y range.
xlogrange: boolean (default = False)
whether the x-axis should have a log scale
ylogrange: boolean (default = False)
whether the y-axis should have a log scale
nlevels: int (default = 20)
number of levels to use for the contours
logscale: boolean (default = True)
whether to use log or linear spaced contours
colorbar: boolean (default = False)
draw a colorbar
scalemin: float (default = arr.min())
minimum value to use for the color scale
scalemax: float (default = arr.max())
maximum value to use for the color scale
pynbody.plot.generic.fourier_map(sim, nbins=100, nmin=1000, nphi=100, mmin=1, mmax=7, rmax=10, levels=[0.01, 0.05, 0.1, 0.2], subplot=None, ret=False, **kwargs)[source]

Plot an overdensity map generated from a Fourier expansion of the particle distribution. A Profile() is made and passed to inv_fourier() to obtain an overdensity map. The map is plotted using the usual matplotlib.contour.

Input:

sim : a SimSnap() object

Optional Keywords:

nbins (100) : number of radial bins to use for the profile

nmin (1000) : minimum number of particles required per bin

nphi (100) : number of azimuthal bins to use for the map

mmin (1) : lowest multiplicity Fourier component

mmax (7) : highest multiplicity Fourier component

rmax (10) : maximum radius to use when generating the profile

levels [0.01,0.05,0.1,0.2] : tuple of levels for plotting contours

subplot (None) : Axes object on which to plot the contours

pynbody.plot.generic.prob_plot(x, y, weight, nbins=(100, 100), extent=None, axes=None, **kwargs)[source]

Make a plot of the probability of y given x, i.e. p(y|x). The values are normalized such that the integral along each column is one.

Input:

x: primary binning axis

y: secondary binning axis

weight: weights array

nbins: tuple of length 2 specifying the number of bins in each direction

extent: tuple of length 4 speciphysical extent of the axes
(xmin,xmax,ymin,ymax)

Optional Keywords:

all optional keywords are passed on to the imshow() command

metals

pynbody.plot.metals.mdf(sim, filename=None, clear=True, range=[-5, 0.3], axes=False, **kwargs)[source]

Metallicity Distribution Function

Plots a metallicity distribution function to the best of matplotlib’s abilities. Unfortunately, the “normed” keyword is buggy and does not return a PDF. The “density” keyword should, but it not yet supported in many versions of numpy.

Usage:

>>> import pynbody.plot as pp
>>> pp.mdf(s,linestyle='dashed',color='k')
pynbody.plot.metals.ofefeh(sim, fxn=<function gauss_kde at 0x108ea5b18>, filename=None, **kwargs)[source]

Use hist2d() to make a [O/Fe] vs. [Fe/H] plot

Input:

sim: snapshot to pull data from

Optional Keywords:

fxn: a function with the same signature as functions
hist2d() and gauss_kde() in pynbody.plot.generic default: pynbody.plot.generic.hist2d()

see make_contour_plot() for other plot-related keywords.

profile

pynbody.plot.profile.rotation_curve(sim, center=True, r_units='kpc', v_units='km s^-1', disk_height='100 pc', nbins=50, bin_spacing='equaln', clear=True, quick=False, filename=None, min=False, max=False, yrange=False, legend=False, parts=False, axes=False, **kwargs)[source]

Centre on potential minimum, align so that the disk is in the x-y plane, then use the potential in that plane to generate and plot a rotation curve.

needs documentation/description of the keyword arguments

pynbody.plot.profile.fourier_profile(sim, center=True, disk_height='2 kpc', nbins=50, pretime='2 Gyr', r_units='kpc', bin_spacing='equaln', clear=True, min=False, max=False, filename=None, **kwargs)[source]

Centre on potential minimum, align so that the disk is in the x-y plane, then plot the amplitude of the 2nd fourier mode as a function of radius.

needs description of the keyword arguments

pynbody.plot.profile.density_profile(sim, linestyle=False, center=True, clear=True, fit=False, in_units=None, filename=None, fit_factor=0.02, axes=False, **kwargs)[source]

3d density profile

Options:

filename (None): name of file to which to save output

Usage:

>>> import pynbody.plot as pp
>>> h = s.halos()
>>> pp.density_profile(h[1],linestyle='dashed',color='k')

sph

routines for plotting smoothed quantities

pynbody.plot.sph.fmt(x, pos)[source]

Custom formatter to handle log of values for color bar Not sure why, but the LogFormatterExponenet sometimes adds and “e” to the value... ?? This simply returns log10 of the value.

pynbody.plot.sph.sideon_image(sim, *args, **kwargs)[source]

Rotate the simulation so that the disc of the passed halo is side-on, then make an SPH image by passing the parameters into the function image

For a description of keyword arguments see image().

pynbody.plot.sph.faceon_image(sim, *args, **kwargs)[source]

Rotate the simulation so that the disc of the passed halo is face-on, then make an SPH image by passing the parameters into the function image

For a description of keyword arguments see image().

pynbody.plot.sph.velocity_image(sim, width='10 kpc', vector_color='black', edgecolor='black', quiverkey_bg_color=None, vector_resolution=40, scale=None, mode='quiver', key_x=0.3, key_y=0.9, key_color='white', key_length='100 km s**-1', quiverkey=True, density=1.0, vector_qty='vel', **kwargs)[source]

Make an SPH image of the given simulation with velocity vectors overlaid on top.

For a description of additional keyword arguments see image(), or see the tutorial.

Keyword arguments:

vector_color (black): The color for the velocity vectors

edgecolor (black): edge color used for the lines - using a color
other than black for the vector_color and a black edgecolor can result in poor readability in pdfs

vector_resolution (40): How many vectors in each dimension (default is 40x40)

quiverkey_bg_color (none): The color for the legend (scale) background

scale (None): The length of a vector that would result in a displayed length of the figure width/height.

mode (‘quiver’): make a ‘quiver’ or ‘stream’ plot

key_x (0.3): Display x (width) position for the vector key (quiver mode only)

key_y (0.9): Display y (height) position for the vector key (quiver mode only)

key_color (white): Color for the vector key (quiver mode only)

key_length (100 km/s): Velocity to use for the vector key (quiver mode only)

density (1.0): Density of stream lines (stream mode only)

quiverkey (True): Whether or not to inset the key

vector_qty (‘vel’): The name of the vector field to plot

pynbody.plot.sph.volume(sim, qty='rho', width=None, resolution=200, color=(1.0, 1.0, 1.0), vmin=None, vmax=None, dynamic_range=4.0, log=True, create_figure=True)[source]

Create a volume rendering of the given simulation using mayavi.

Keyword arguments:

qty (rho): The name of the array to interpolate

width (None): The width of the cube to generate, centered on the origin

resolution (200): The number of elements along each side of the cube

color (white): The color of the volume rendering. The value of each voxel
is used to set the opacity.

vmin (None): The value for zero opacity (calculated using dynamic_range if None)

vmax (None): The value for full opacity (calculated from the maximum
value in the region if None)

dynamic_range: The dynamic range to use if vmin and vmax are not specified

log (True): log-scale the image before passing to mayavi

create_figure (True): create a new mayavi figure before rendering

pynbody.plot.sph.contour(*args, **kwargs)[source]

Make an SPH image of the given simulation and render it as contours. nlevels and levels are passed to pyplot’s contour command.

Other arguments are as for image.

pynbody.plot.sph.image(sim, qty='rho', width='10 kpc', resolution=500, units=None, log=True, vmin=None, vmax=None, av_z=False, filename=None, z_camera=None, clear=True, cmap=None, title=None, qtytitle=None, show_cbar=True, subplot=False, noplot=False, ret_im=False, fill_nan=True, fill_val=0.0, linthresh=None, **kwargs)[source]

Make an SPH image of the given simulation.

Keyword arguments:

qty (rho): The name of the array to interpolate

width (10 kpc): The overall width and height of the plot. If
width is a float or an int, then it is assumed to be in units of sim['pos']. It can also be passed in as a string indicating the units, i.e. ‘10 kpc’, in which case it is converted to units of sim['pos'].

resolution (500): The number of pixels wide and tall

units (None): The units of the output

av_z (False): If True, the requested quantity is averaged down
the line of sight (default False: image is generated in the thin plane z=0, unless output units imply an integral down the line of sight). If a string, the requested quantity is averaged down the line of sight weighted by the av_z array (e.g. use ‘rho’ for density-weighted quantity; the default results when av_z=True are volume-weighted).
z_camera (None): If set, a perspective image is rendered. See
pynbody.sph.image() for more details.

filename (None): if set, the image will be saved in a file

clear (True): whether to call clf() on the axes first

cmap (None): user-supplied colormap instance

title (None): plot title

qtytitle (None): colorbar quantity title

show_cbar (True): whether to plot the colorbar

subplot (False): the user can supply a AxesSubPlot instance on which the image will be shown

noplot (False): do not display the image, just return the image array

ret_im (False): return the image instance returned by imshow

num_threads (None) : if set, specify the number of threads for the multi-threaded routines; otherwise the pynbody.config default is used

fill_nan (True): if any of the image values are NaN, replace with fill_val

fill_val (0.0): the fill value to use when replacing NaNs

linthresh (None): if the image has negative and positive values
and a log scaling is requested, the part between -linthresh and linthresh is shown on a linear scale to avoid divergence at 0

stars

pynbody.plot.stars.render(sim, filename=None, r_band='i', g_band='v', b_band='u', r_scale=0.5, g_scale=1.0, b_scale=1.0, dynamic_range=2.0, mag_range=None, width=50, starsize=None, plot=True, axes=None, ret_im=False, clear=True, ret_range=False)[source]

Make a 3-color image of stars.

The colors are based on magnitudes found using stellar Marigo stellar population code. However there is no radiative transfer to account for dust.

Returns: If ret_im=True, an NxNx3 array representing an RGB image

Optional keyword arguments:

filename: string (default: None)
Filename to be written to (if a filename is specified)
r_band: string (default: ‘i’)
Determines which Johnston filter will go into the image red channel
g_band: string (default: ‘v’)
Determines which Johnston filter will go into the image green channel
b_band: string (default: ‘b’)
Determines which Johnston filter will go into the image blue channel
r_scale: float (default: 0.5)
The scaling of the red channel before channels are combined
g_scale: float (default: 1.0)
The scaling of the green channel before channels are combined
b_scale: float (default: 1.0)
The scaling of the blue channel before channels are combined
width: float in kpc (default:50)
Sets the size of the image field in kpc
starsize: float in kpc (default: None)
If not None, sets the maximum size of stars in the image
ret_im: bool (default: False)
if True, the NxNx3 image array is returned
ret_range: bool (default: False)
if True, the range of the image in mag arcsec^-2 is returned.
plot: bool (default: True)
if True, the image is plotted
axes: matplotlib axes object (deault: None)
if not None, the axes object to plot to
dynamic_range: float (default: 2.0)
The number of dex in luminosity over which the image brightness ranges
mag_range: float, float (default: None)
If provided, the brightest and faintest surface brightnesses in the range, in mag arcsec^-2. Takes precedence over dynamic_range.
pynbody.plot.stars.sfh(sim, filename=None, massform=True, clear=False, legend=False, subplot=False, trange=False, bins=100, **kwargs)[source]

star formation history

Optional keyword arguments:

trange: list, array, or tuple
size(t_range) must be 2. Specifies the time range.
bins: int
number of bins to use for the SFH
massform: bool
decides whether to use original star mass (massform) or final star mass
subplot: subplot object
where to plot SFH
legend: boolean
whether to draw a legend or not
clear: boolean
if False (default), plot on the current axes. Otherwise, clear the figure first.

By default, sfh will use the formation mass of the star. In tipsy, this will be taken from the starlog file. Set massform=False if you want the final (observed) star formation history

Usage:

>>> import pynbody.plot as pp
>>> pp.sfh(s,linestyle='dashed',color='k')
pynbody.plot.stars.schmidtlaw(sim, center=True, filename=None, pretime='50 Myr', diskheight='3 kpc', rmax='20 kpc', compare=True, radial=True, clear=True, legend=True, bins=10, **kwargs)[source]

Schmidt Law

Plots star formation surface density vs. gas surface density including the observed relationships. Currently, only plots densities found in radial annuli.

Usage:

>>> import pynbody.plot as pp
>>> pp.schmidtlaw(h[1])

Optional keyword arguments:

center: bool
center and align the input simulation faceon.
filename: string
Name of output file

pretime (default=‘50 Myr’): age of stars to consider for SFR

diskheight (default=‘3 kpc’): height of gas and stars above
and below disk considered for SF and gas densities.

rmax (default=‘20 kpc’): radius of disk considered

compare (default=True): whether to include Kennicutt (1998) and
Bigiel+ (2008) for comparison

radial (default=True): should bins be annuli or a rectangular grid?

bins (default=10): How many radial bins should there be?

legend: boolean
whether to draw a legend or not
pynbody.plot.stars.oneschmidtlawpoint(sim, center=True, pretime='50 Myr', diskheight='3 kpc', rmax='20 kpc', **kwargs)[source]

One Schmidt Law Point

Determines values for star formation surface density and gas surface density for the entire galaxy based on the half mass cold gas radius.

Usage: import pynbody.plot as pp pp.oneschmidtlawpoint(h[1])

pretime (default=‘50 Myr’): age of stars to consider for SFR

diskheight (default=‘3 kpc’): height of gas and stars above
and below disk considered for SF and gas densities.

rmax (default=‘20 kpc’): radius of disk considered

pynbody.plot.stars.satlf(sim, band='v', filename=None, MWcompare=True, Trentham=True, clear=True, legend=True, label='Simulation', **kwargs)[source]

satellite luminosity function

Options:

band (‘v’): which Johnson band to use. available filters: U, B, V, R, I, J, H, K

filename (None): name of file to which to save output

MWcompare (True): whether to plot comparison lines to MW

Trentham (True): whether to plot comparison lines to Trentham +
Tully (2009) combined with Koposov et al (2007)

By default, satlf will use the formation mass of the star. In tipsy, this will be taken from the starlog file.

Usage:

>>> import pynbody.plot as pp
>>> h = s.halos()
>>> pp.satlf(h[1],linestyle='dashed',color='k')
pynbody.plot.stars.sbprofile(sim, band='v', diskheight='3 kpc', rmax='20 kpc', binning='equaln', center=True, clear=True, filename=None, axes=False, fit_exp=False, print_ylabel=True, fit_sersic=False, **kwargs)[source]

surface brightness profile

Usage:

>>> import pynbody.plot as pp
>>> h = s.halos()
>>> pp.sbprofile(h[1],exp_fit=3,linestyle='dashed',color='k')

Options:

band (‘v’): which Johnson band to use. available filters: U, B,
V, R, I, J, H, K

*fit_exp*(False): Fits straight exponential line outside radius specified.

*fit_sersic*(False): Fits Sersic profile outside radius specified.

diskheight(‘3 kpc’) rmax(‘20 kpc’): Size of disk to be profiled

binning(‘equaln’): How show bin sizes be determined? based on
pynbody.analysis.profile

center(True): Automatically align face on and center?

axes(False): In which axes (subplot) should it be plotted?

filename (None): name of file to which to save output

needs a description of all keywords

By default, sbprof will use the formation mass of the star. In tipsy, this will be taken from the starlog file.

pynbody.plot.stars.moster(xmasses, z)[source]

Based on Moster+ (2013) return what stellar mass corresponds to the halo mass passed in.

Usage

>>> from pynbody.plot.stars import moster
>>> xmasses = np.logspace(np.log10(min(totmasshalos)),1+np.log10(max(totmasshalos)),20)
>>> ystarmasses, errors = moster(xmasses,halo_catalog._halos[1].properties['z'])
>>> plt.fill_between(xmasses,np.array(ystarmasses)/np.array(errors),
                  y2=np.array(ystarmasses)*np.array(errors),
                  facecolor='#BBBBBB',color='#BBBBBB')
pynbody.plot.stars.behroozi(xmasses, z)[source]

Based on Behroozi+ (2013) return what stellar mass corresponds to the halo mass passed in.

Usage

>>> from pynbody.plot.stars import moster
>>> xmasses = np.logspace(np.log10(min(totmasshalos)),1+np.log10(max(totmasshalos)),20)
>>> ystarmasses, errors = moster(xmasses,halo_catalog._halos[1].properties['z'])
>>> plt.fill_between(xmasses,np.array(ystarmasses)/np.array(errors),
                  y2=np.array(ystarmasses)*np.array(errors),
                  facecolor='#BBBBBB',color='#BBBBBB')
pynbody.plot.stars.subfindguo(halo_catalog, clear=False, compare=True, baryfrac=False, filename=False, **kwargs)[source]

Stellar Mass vs. Halo Mass

Takes a halo catalogue and plots the member stellar masses as a function of halo mass.

Usage:

>>> import pynbody.plot as pp
>>> h = s.halos()
>>> pp.guo(h,marker='+',markerfacecolor='k')

Options:

compare (True): Should comparison line be plotted?
If compare = ‘guo’, Guo+ (2010) plotted instead of Behroozi+ (2013)

baryfrac (False): Should line be drawn for cosmic baryon fraction?

filename (None): name of file to which to save output

pynbody.plot.stars.guo(halo_catalog, clear=False, compare=True, baryfrac=False, filename=False, **kwargs)[source]

Stellar Mass vs. Halo Mass

Takes a halo catalogue and plots the member stellar masses as a function of halo mass.

Usage:

>>> import pynbody.plot as pp
>>> h = s.halos()
>>> pp.guo(h,marker='+',markerfacecolor='k')

Options:

compare (True): Should comparison line be plotted?
If compare = ‘guo’, Guo+ (2010) plotted instead of Behroozi+ (2013)

baryfrac (False): Should line be drawn for cosmic baryon fraction?

filename (None): name of file to which to save output

util

Utility functions for the plotting module

pynbody.plot.util.fast_kde(x, y, kern_nx=None, kern_ny=None, gridsize=(100, 100), extents=None, nocorrelation=False, weights=None, norm=False, **kwargs)[source]

A faster gaussian kernel density estimate (KDE). Intended for computing the KDE on a regular grid (different use case than scipy’s original scipy.stats.kde.gaussian_kde()).

Author: Joe Kington License: MIT License <http://www.opensource.org/licenses/mit-license.php>

Performs a gaussian kernel density estimate over a regular grid using a convolution of the gaussian kernel with a 2D histogram of the data.

This function is typically several orders of magnitude faster than scipy.stats.kde.gaussian_kde for large (>1e7) numbers of points and produces an essentially identical result.

Input:

x: array
The x-coords of the input data points
y: array
The y-coords of the input data points
kern_nx: float
size (in units of x) of the kernel
kern_ny: float
size (in units of y) of the kernel
gridsize: (Nx , Ny) tuple (default: 200x200)
Size of the output grid
extents: (default: extent of input data) A (xmin, xmax, ymin, ymax)
tuple of the extents of output grid
nocorrelation: (default: False) If True, the correlation between the
x and y coords will be ignored when preforming the KDE.
weights: (default: None) An array of the same shape as x & y that
weighs each sample (x_i, y_i) by each value in weights (w_i). Defaults to an array of ones the same size as x & y.
norm: boolean (default: False)
If False, the output is only corrected for the kernel. If True, the result is normalized such that the integral over the area yields 1.
Output:
A gridded 2D kernel density estimate of the input points.
pynbody.plot.util.inv_fourier(p, nmin=1000, mmin=1, mmax=7, nphi=100)[source]

Invert a profile with fourier coefficients to yield an overdensity map.

Inputs:

p : a Profile() object

Optional Keywords:

nmin (1000) : minimum number of particles required per bin

mmin (1) : lowest multiplicity Fourier component

mmax (7) : highest multiplicity Fourier component

nphi (100) : number of azimuthal bins to use for the map