Generic Convenience Modules

bridge

The bridge module has tools for connecting different outputs. For instance, it’s possible to take a subview (e.g. a halo) from one snapshot and ‘push’ it into the other. This is especially useful if the two snapshots are different time outputs of the same simulation.

Once connected, bridge called on a specific subset of particles in output1 will trace these particles back (or forward) to the output2, enabling observing a change in their properties, such as position, temperature, etc.

For a tutorial on how to use the bridge module to trace the particles in your simulation, see the bridge tutorial.

class pynbody.bridge.Bridge(start, end)[source]

Generic Bridge class

Methods

__call__(s)

Given a subview of either the start or end point of the bridge, generate the corresponding subview of the connected snapshot

catalog_transfer_matrix([min_index, …])

Return a max_index x max_index matrix with the number of particles transferred from the row group in groups_1 to the column group in groups_2.

fuzzy_match_catalog([min_index, max_index, …])

fuzzy_match_catalog returns, for each halo in groups_1, a list of possible identifications in groups_2, along with the fraction of particles in common between the two.

is_same(i2)

Returns true if the particle i1 in the start point is the same as the particle i2 in the end point.

match_catalog([min_index, max_index, …])

Given a Halos object groups_1, a Halos object groups_2 and a Bridge object connecting the two parent simulations, this identifies the most likely ID’s in groups_2 of the objects specified in groups_1.

is_same(i2)[source]

Returns true if the particle i1 in the start point is the same as the particle i2 in the end point.

match_catalog(min_index=1, max_index=30, threshold=0.5, groups_1=None, groups_2=None, use_family=None)[source]

Given a Halos object groups_1, a Halos object groups_2 and a Bridge object connecting the two parent simulations, this identifies the most likely ID’s in groups_2 of the objects specified in groups_1.

If groups_1 and groups_2 are not specified, they are automatically obtained using the SimSnap.halos method.

Parameters min_index and max_index are the minimum and maximum halo numbers to be matched (in both ends of the bridge). If max_index is too large, the catalogue matching can take prohibitively long (it scales as max_index^2).

This routine currently uses particle number as a proxy for mass, so that the main simulation data does not need to be loaded.

If b links snapshot f1 (high redshift) to f2 (low redshift) and we set

cat = b.match_catalog()

then cat is now a numpy index array such that f1.halos()[i] is the major progenitor for f2.halos()[cat[i]], assuming cat[i] is positive.

cat[0:min_index+1] is set to -2. Halos which cannot be matched because they have too few particles in common give the result -1. This is determined by the given threshold fraction of particles in common (by default, 50%).

If use_family is specified, only particles from that family are cross-matched. This can be useful e.g. if matching between two different simulations where the relationship between DM particles is known, but perhaps the relationship between star particles is random.

fuzzy_match_catalog(min_index=1, max_index=30, threshold=0.01, groups_1=None, groups_2=None, use_family=None, only_family=None)[source]

fuzzy_match_catalog returns, for each halo in groups_1, a list of possible identifications in groups_2, along with the fraction of particles in common between the two.

Normally, match_catalog is simpler to use, but this routine offers greater flexibility for advanced users. The first entry for each halo corresponds to the output from match_catalog.

If no identification is found, the entry is the empty list [].

catalog_transfer_matrix(min_index=1, max_index=30, groups_1=None, groups_2=None, use_family=None, only_family=None)[source]

Return a max_index x max_index matrix with the number of particles transferred from the row group in groups_1 to the column group in groups_2.

Normally, match_catalog (or fuzzy_match_catalog) are easier to use, but this routine provides the maximal information.

class pynbody.bridge.OrderBridge(start, end, order_array='iord', monotonic=True, allow_family_change=False)[source]

An OrderBridge uses integer arrays in two simulations (start,end) where particles i_start and i_end are defined to be the same if and only if start[order_array][i_start] == start[order_array][i_end].

If monotonic is True, order_array must be monotonically increasing in both ends of the bridge (and this is not checked for you). If monotonic is False, the bridging is slower but this is the failsafe option.

Methods

__call__(s)

Given a subview of either the start or end point of the bridge, generate the corresponding subview of the connected snapshot

catalog_transfer_matrix([min_index, …])

Return a max_index x max_index matrix with the number of particles transferred from the row group in groups_1 to the column group in groups_2.

fuzzy_match_catalog([min_index, max_index, …])

fuzzy_match_catalog returns, for each halo in groups_1, a list of possible identifications in groups_2, along with the fraction of particles in common between the two.

is_same(i1, i2)

Returns true if the particle i1 in the start point is the same as the particle i2 in the end point.

match_catalog([min_index, max_index, …])

Given a Halos object groups_1, a Halos object groups_2 and a Bridge object connecting the two parent simulations, this identifies the most likely ID’s in groups_2 of the objects specified in groups_1.

is_same(i1, i2)[source]

Returns true if the particle i1 in the start point is the same as the particle i2 in the end point.

pynbody.bridge.bridge_factory(a, b)[source]

Create a bridge connecting the two specified snapshots. For more information see bridge-tutorial.

filt

Defines and implements ‘filters’ which allow abstract subsets of data to be specified.

See the filter tutorial for some sample usage.

class pynbody.filt.Sphere(radius, cen=0, 0, 0)[source]

Return particles that are within radius of the point cen.

Methods

__call__(sim)

Call self as a function.

where

class pynbody.filt.Cuboid(x1, y1=None, z1=None, x2=None, y2=None, z2=None)[source]

Create a cube with specified edge coordinates. If any of the cube coordinates x1, y1, z1, x2, y2, z2 are not specified they are determined as y1=x1; z1=x1; x2=-x1; y2=-y1; z2=-z1.

Methods

__call__(sim)

Call self as a function.

where

class pynbody.filt.Disc(radius, height, cen=0, 0, 0)[source]

Return particles that are within a disc of extent radius and thickness height centered on cen.

Methods

__call__(sim)

Call self as a function.

where

class pynbody.filt.BandPass(prop, min, max)[source]

Return particles whose property prop is within min and max, which can be specified as unit strings.

Methods

__call__(sim)

Call self as a function.

where

class pynbody.filt.HighPass(prop, min)[source]

Return particles whose property prop exceeds min, which can be specified as a unit string.

Methods

__call__(sim)

Call self as a function.

where

class pynbody.filt.LowPass(prop, max)[source]

Return particles whose property prop is less than max, which can be specified as a unit string.

Methods

__call__(sim)

Call self as a function.

where

pynbody.filt.Annulus(r1, r2, cen=0, 0, 0)[source]

Convenience function that returns a filter which selects particles in between two spheres specified by radii r1 and r2 centered on cen.

pynbody.filt.SolarNeighborhood(r1=Unit('5.00e+00 kpc'), r2=Unit('1.00e+01 kpc'), height=Unit('2.00e+00 kpc'), cen=0, 0, 0)[source]

Convenience function that returns a filter which selects particles in a disc between radii r1 and r2 and thickness height.

halo

Implements halo catalogue functions. If you have a supported halo catalogue on disk or a halo finder installed and correctly configured, you can access a halo catalogue through f.halos() where f is a SimSnap.

See the halo tutorial for some examples.

class pynbody.halo.Halo(halo_id, halo_catalogue, *args, **kwa)[source]

Generic class representing a halo.

Attributes
ancestor

The original SimSnap from which this view is derived (potentially self)

filename

Methods

all_keys()

Returns a list of all arrays that can be either lazy-evaluated or lazy loaded from an auxiliary file.

bridge(other)

Tries to construct a bridge function between this SimSnap and another one.

conversion_context()

Return a dictionary containing a (scalefactor) and h (Hubble constant in canonical units) for this snapshot, ready for passing into unit conversion functions.

derivable_keys()

Returns a list of arrays which can be lazy-evaluated.

derived_quantity(fn)

families()

Return the particle families which have representitives in this SimSnap.

family_keys([fam])

Return list of arrays which are not accessible from this view, but can be accessed from family-specific sub-views.

get(key[, alternative])

Standard python get method, returns self[key] if key in self else alternative

get_index_list(relative_to[, of_particles])

Get a list specifying the index of the particles in this view relative to the ancestor relative_to, such that relative_to[get_index_list(relative_to)]==self.

halos(*args, **kwargs)

Tries to instantiate a halo catalogue object for the given snapshot, using the first available method (as defined in the configuration files).

has_family_key(name)

Returns True if the array name is accessible (in memory) for at least one family

has_key(name)

Returns True if the array name is accessible (in memory)

infer_original_units(*args)

Return the units on disk for a quantity with the specified dimensions

intersect(other[, op])

Returns the set intersection of this simulation view with another view of the same simulation

is_ancestor(other)

Returns true if other is a subview of self

is_derived_array(v[, fam])

Returns True if the array or family array of given name is auto-derived (and therefore read-only).

is_descendant(other)

Returns true if self is a subview of other

is_subhalo(otherhalo)

Convenience function that calls the corresponding function in a halo catalogue.

items()

Returns a list of tuples describing the array names and their contents in memory

keys()

Return the directly accessible array names (in memory)

load_copy()

Tries to load a copy of this snapshot, using partial loading to select only a subset of particles corresponding to a given SubSnap

loadable_keys([fam])

Returns a list of arrays which can be lazy-loaded from an auxiliary file.

mean_by_mass(name)

Calculate the mean by mass of the specified array.

original_units()

Converts all arrays’units to be consistent with the units of the original file.

physical_units(*args, **kwargs)

Converts all array’s units to be consistent with the distance, velocity, mass basis units specified.

rotate_x(angle)

Rotates the snapshot about the current x-axis by ‘angle’ degrees.

rotate_y(angle)

Rotates the snapshot about the current y-axis by ‘angle’ degrees.

rotate_z(angle)

Rotates the snapshot about the current z-axis by ‘angle’ degrees.

set_units_system([velocity, distance, mass, …])

Set the unit system for the snapshot by specifying any or all of velocity, distance, mass and temperature units.

setdiff(other)

Returns the set difference of this simulation view with another view of the same simulation

transform(matrix)

union(other)

Returns the set union of this simulation view with another view of the same simulation

unlink_array(name)

If the named array is auto-derived, this destroys the link so that the array becomes editable but no longer auto-updates.

values()

Returns a list of the actual arrays in memory

wrap([boxsize, convention])

Wraps the positions of the particles in the box to lie between [-boxsize/2, boxsize/2].

write([fmt, filename])

write_array(array_name[, fam])

Write out the array with the specified name.

decorator

iteritems

iterkeys

itervalues

stable_derived_quantity

is_subhalo(otherhalo)[source]

Convenience function that calls the corresponding function in a halo catalogue.

class pynbody.halo.HaloCatalogue(sim)[source]

Generic halo catalogue object.

Attributes
base

Methods

get_group_array()

Return an array with an integer for each particle in the simulation indicating which halo that particle is associated with.

is_subhalo(childid, parentid)

Checks whether the specified ‘childid’ halo is a subhalo of ‘parentid’ halo.

calc_item

contains

is_subhalo(childid, parentid)[source]

Checks whether the specified ‘childid’ halo is a subhalo of ‘parentid’ halo.

get_group_array()[source]

Return an array with an integer for each particle in the simulation indicating which halo that particle is associated with. If there are multiple levels (i.e. subhalos), the number returned corresponds to the lowest level, i.e. the smallest subhalo.

class pynbody.halo.GrpCatalogue(sim, array='grp', ignore=None, **kwargs)[source]

A generic catalogue using a .grp file to specify which particles belong to which group.

Attributes
base

Methods

get_group_array([family])

Return an array with an integer for each particle in the simulation indicating which halo that particle is associated with.

is_subhalo(childid, parentid)

Checks whether the specified ‘childid’ halo is a subhalo of ‘parentid’ halo.

load_copy(i)

Load the a fresh SimSnap with only the particle in halo i

precalculate()

Speed up future operations by precalculating the indices for all halos in one operation.

calc_item

contains

precalculate()[source]

Speed up future operations by precalculating the indices for all halos in one operation. This is slow compared to getting a single halo, however.

get_group_array(family=None)[source]

Return an array with an integer for each particle in the simulation indicating which halo that particle is associated with. If there are multiple levels (i.e. subhalos), the number returned corresponds to the lowest level, i.e. the smallest subhalo.

load_copy(i)[source]

Load the a fresh SimSnap with only the particle in halo i

class pynbody.halo.AmigaGrpCatalogue(sim, arr_name='amiga.grp', **kwargs)[source]
Attributes
base

Methods

get_group_array([family])

Return an array with an integer for each particle in the simulation indicating which halo that particle is associated with.

is_subhalo(childid, parentid)

Checks whether the specified ‘childid’ halo is a subhalo of ‘parentid’ halo.

load_copy(i)

Load the a fresh SimSnap with only the particle in halo i

precalculate()

Speed up future operations by precalculating the indices for all halos in one operation.

calc_item

contains

sph

pynbody SPH rendering module.

This module encompasses Kernel objects, which return C fragments from which a final C code to perform the rendering is derived.

For most users, the function of interest will be render_image().

pynbody.sph.render_spherical_image(snap, qty='rho', nside=8, distance=10.0, kernel=<pynbody.sph.Kernel object>, kstep=0.5, denoise=None, out_units=None, threaded=None)[source]

Render an SPH image on a spherical surface. Requires healpy libraries.

Keyword arguments:

qty (‘rho’): The name of the simulation array to render

nside (8): The healpix nside resolution to use (must be power of 2)

distance (10.0): The distance of the shell (for 3D kernels) or maximum distance

of the skewers (2D kernels)

kernel: The Kernel object to use (defaults to 3D spline kernel)

kstep (0.5): The sampling distance when projecting onto the spherical surface in units of the

smoothing length

denoise (False): if True, divide through by an estimate of the discreteness noise.

The returned image is then not strictly an SPH estimate, but this option can be useful to reduce noise.

threaded: if False, render on a single core. Otherwise, the number of threads to use.

Defaults to a value specified in your configuration files. Currently multi-threaded rendering is slower than single-threaded because healpy does not release the gil.

pynbody.sph.render_image(snap, qty='rho', x2=100, nx=500, y2=None, ny=None, x1=None, y1=None, z_plane=0.0, out_units=None, xy_units=None, kernel=<pynbody.sph.Kernel object>, z_camera=None, smooth='smooth', smooth_in_pixels=False, force_quiet=False, approximate_fast=True, threaded=None, denoise=None)[source]

Render an SPH image using a typical (mass/rho)-weighted ‘scatter’ scheme.

Keyword arguments:

qty (‘rho’): The name of the array within the simulation to render

x2 (100.0): The x-coordinate of the right edge of the image

nx (500): The number of pixels wide to make the image

y2: The y-coordinate of the upper edge of the image (default x2,

or if ny is specified, x2*ny/nx)

ny (nx): The number of pixels tall to make the image

x1 (-x2): The x-coordinate of the left edge of the image

y1 (-y2): The y-coordinate of the lower edge of the image

z_plane (0.0): The z-coordinate of the plane of the image

out_units (no conversion): The units to convert the output image into

xy_units: The units for the x and y axes

kernel: The Kernel object to use (default Kernel(), a 3D spline kernel)

z_camera: If this is set, a perspective image is rendered,

assuming the kernel is suitable (i.e. is a projecting kernel). The camera is at the specified z coordinate looking towards -ve z, and each pixel represents a line-of-sight radially outwards from the camera. The width then specifies the width of the image in the z=0 plane. Particles too close to the camera are also excluded.

smooth: The name of the array which contains the smoothing lengths

(default ‘smooth’)

smooth_in_pixels: If True, the smoothing array contains the smoothing

length in image pixels, rather than in real distance units (default False)

approximate_fast: if True, render high smoothing length particles at

progressively lower resolution, resample and sum

denoise: if True, divide through by an estimate of the discreteness noise.

The returned image is then not strictly an SPH estimate, but this option can be useful to reduce noise especially when rendering AMR grids which often introduce problematic edge effects.

verbose: if True, all text output suppressed

threaded: if False (or None), render on a single core. Otherwise,

the number of threads to use (defaults to a value specified in your configuration files).

pynbody.sph.to_3d_grid(snap, qty='rho', nx=None, ny=None, nz=None, x2=None, out_units=None, xy_units=None, kernel=<pynbody.sph.Kernel object>, smooth='smooth', approximate_fast=True, threaded=None, snap_slice=None, denoise=None)[source]

Project SPH onto a grid using a typical (mass/rho)-weighted ‘scatter’ scheme.

Keyword arguments:

qty (‘rho’): The name of the array within the simulation to render

nx (x2-x1 / soft): The number of pixels wide to make the grid

ny (nx): The number of pixels tall to make the grid

nz (nx): The number of pixels deep to make the grid

out_units (no conversion): The units to convert the output grid into

xy_units: The units for the x and y axes

kernel: The Kernel object to use (default Kernel(), a 3D spline kernel)

smooth: The name of the array which contains the smoothing lengths

(default ‘smooth’)

denoise: if True, divide through by an estimate of the discreteness noise.

The returned image is then not strictly an SPH estimate, but this option can be useful to reduce noise especially when rendering AMR grids which often introduce problematic edge effects.

pynbody.sph.spectra(snap, qty='rho', x1=0.0, y1=0.0, v2=400, nvel=200, v1=None, element='H', ion='I', xy_units=Unit('kpc'), vel_units=Unit('km s**-1'), smooth='smooth', __threaded=False)[source]

Render an SPH spectrum using a (mass/rho)-weighted ‘scatter’ scheme of all the particles that have a smoothing length within 2 h_sm of the position.

Keyword arguments:

qty (‘rho’): The name of the array within the simulation to render

x1 (0.0): The x-coordinate of the line of sight.

y1 (0.0): The y-coordinate of the line of sight.

v1 (-400.0): The minimum velocity of the spectrum

v2 (400.0): The maximum velocity of the spectrum

nvel (500): The number of resolution elements in spectrum

xy_units (‘kpc’): The units for the x and y axes

smooth: The name of the array which contains the smoothing lengths

(default ‘smooth’)