# Author: Volker Hoffmann <volmker@cheleb.net>, <volker@physik.uzh.ch>
# Date: 02 December 2014
%matplotlib inline
import matplotlib as mpl; mpl.rcParams['savefig.dpi'] = 144
import numpy as np
import matplotlib.pyplot as plt
Using Pynbody & Pkdgrav2 Data¶
In this notebook, we will use Pynbody [1] to load particle data from Pkdgrav2 [2]. Pkdgrav2 is a dark-matter only simulation, so we will only be using a small subset of Pynbody can do. Using different backends, we can, for example, load SPH particle data and use Pynbody to construct density, velocity, and pressure fields. We will expore this in a forthcoming article. Note that this notebook is somewhat compressed, and we refer to the official Pynbody introduction [3] for more details.
import pynbody
Load File Data¶
To load Pynbody simulation data, all we need to do is point to the simulation output. Note that Pynbody lazy-loads, i.e. it doesn't actually load the data unless they are specifically requested.
# Load Data, Show Some Stuff
s = pynbody.load('/zbox/data/dpotter/sr8.00100')
print s
This simulation has only dark matter particles. Let's first check how many, and then check some properties, such as units.
print len(s)
print len(s.dark)
print s.properties
We mentioned that Pynbody lazy-loads. So if we check what data are available to us with keys(), we will get an empty list because nothing has actually been loaded. We can use loadable_keys() to check what data we can load. We can also call all_keys() to check all the stuff that could be avaible through some post-processing routines.
print s.keys()
print s.loadable_keys()
print s.all_keys()
Directly Accessing Particle Data¶
For example, to get a arrays with the mass, positions, and velcoity components of all dark matter particles, we can do the following. We can also access the gravitational softening of a particular particle ("eps") or the local gravitational potential ("phi").
print s.dark["mass"]
print s.dark["pos"]
print s.dark["vel"]
Custom Plots¶
With have access to the raw data, we can build any analysis and plots we want. For example, a histogram in mass could be made like follows. Note that this isn't a particularly good example though.
fig, ax = plt.subplots(1,1)
# Pretty Custom Color
teal = ( 119./255., \
151./255., \
161./255. )
# ax.hist(np.log10(i[a>0]), bins=128, log=True, color=c_amazing, histtype="stepfilled", alpha=0.5)
ax.hist(np.log10(s.dark["mass"]), bins=8, log=True, \
color=teal, histtype="bar", \
alpha=0.5)
ax.set_title("Mass Function")
ax.set_xlabel("Log10(Mass)")
ax.set_ylabel("dN/dM")
pass
Pynbody Plots¶
We can also use Pynbody's extensive plotting and analysis routines to make much prettier plots.
im = pynbody.plot.image(s.dark, width='1 kpc', \
units="Msol kpc^-2", cmap=mpl.cm.bone)
