Slides adapted from Britton Smith's Python in Astronomy 2015 talk¶
Data on disk has no inherent physical meaning¶
yt lets you think about the data using a physically motivated interface¶
yt lets you think about the data using a physically motivated interface¶
In [1]:
import yt
#
# download from:
# https://hub.yt/data/goldbaum2016a/feedback_20pc/simulation_outputs/DD0600.tar.gz
# (3.5 GB)
#
ds = yt.load('DD0600/DD0600')
Allowing you to forget about what your data looks like as a file format¶
And select only the data you want to select¶
In [2]:
from yt.units import kpc
sp = ds.sphere(ds.domain_center, 20*kpc)
yt natively deals with multiresolution data¶
Data from each selected cell is returned as a flat, 1D array, with unit metadata attached¶
In [3]:
print(sp['density'])
Data objects can be queried for many fields¶
In [4]:
print(sp['temperature'])
Data objects can be queried for many fields¶
In [5]:
from pprint import pprint
# fields that are defined in the on-disk data file
pprint(ds.field_list)
Data objects can be queried for many fields¶
In [6]:
# fields that yt can calculate given the fields in ds.field_list
# (only showing gas fields to fit on one screen)
pprint([f for f in ds.derived_field_list if f[0]=='gas'])
Spatial information is not lost¶
In [7]:
sp['x']
Out[7]:
In [8]:
sp['x'].to('kpc')
Out[8]:
Spatial information is not lost¶
In [9]:
import numpy as np
np.unique(sp['dx']).to('kpc')
Out[9]:
Data objects return fields as unit-aware numpy arrays¶
In [10]:
from yt.units import kboltz, mh
pressure = sp['density']*sp['temperature']*kboltz/(0.6*mh)
pressure.in_units('Pa')
Out[10]:
In [11]:
pressure.in_units('dyne/cm**2')
Out[11]:
Derived quantities turn fields into single values¶
In [12]:
sp['cell_mass']
Out[12]:
Derived quantities turn fields into single values¶
In [13]:
sp.quantities.total_quantity('cell_mass')
Out[13]:
In [14]:
sp.min('temperature')
Out[14]:
In [15]:
sp.max('temperature')
Out[15]:
In [16]:
sp.mean('temperature', weight='cell_mass')
Out[16]:
In [17]:
sp.mean('temperature', weight='cell_volume')
Out[17]:
In [18]:
from IPython.display import IFrame
IFrame("http://yt-project.org/docs/dev/analyzing/objects.html#available-derived-quantities", width=700, height=600)
Out[18]:
Create new fields by defining a python function¶
In [19]:
def my_entropy(field, data):
return (kboltz * data['temperature'] / data['number_density']**(2/3))
ds.add_field('entropy', function=my_entropy, units="keV*cm**2")
In [20]:
sp['entropy']
Out[20]:
Data Containers¶
In [21]:
ad = ds.all_data()
box = ds.box(ds.domain_left_edge + 500*kpc, ds.domain_right_edge - 500*kpc)
sp = ds.sphere(ds.domain_center, 500*kpc)
disk = ds.disk(ds.domain_center, [0, 0, 1], 500*kpc, 100*kpc)
ray = ds.ray(ds.domain_left_edge, ds.domain_right_edge)
In [22]:
IFrame("http://yt-project.org/doc/analyzing/objects.html#available-objects", width=700, height=600)
Out[22]:
Particle filters¶
In [23]:
from yt.data_objects.particle_filters import add_particle_filter
def young_star(pfilter, data):
filter = data.ds.current_time - data[pfilter.filtered_type, "creation_time"]
filter.convert_to_units('Myr')
return filter < 10
add_particle_filter("young_stars", function=young_star, filtered_type='io',
requires=["creation_time"])
In [24]:
ds.add_particle_filter('young_stars')
ad = ds.all_data()
print(ad['io', 'particle_mass'].shape)
print(ad['young_stars', 'particle_mass'].shape)
This uses ParticlePlot, a new feature in yt 3.3
In [25]:
from yt import ParticlePlot
p = ParticlePlot(ds,
('young_stars', 'particle_position_x'),
('young_stars', 'particle_position_y'),
('young_stars', 'particle_mass'),
width=(10, 'kpc'))
p.set_unit(('young_stars', 'particle_mass'), 'Msun')
Out[25]:
In [26]:
from yt import ProjectionPlot
prj = ProjectionPlot(ds, 2, ('gas', 'density'), width=(15, 'kpc'))
In [27]:
dd = ds.sphere(ds.domain_center, (10, 'kpc'))
prj = ProjectionPlot(ds, 2, 'density', data_source=dd, width=(20, 'kpc'))
prj.set_zlim('density', 1e-4, 1e0)
Out[27]:
In [28]:
ProjectionPlot(ds, 2, 'temperature', weight_field='density', width=(15, 'kpc'))
Out[28]:
In [29]:
import numexpr as ne
from yt.units import G
from numpy import pi
def _total_dynamical_time(field, data):
dens = data['density']
dmd = data['dark_matter_density']
dens.convert_to_cgs()
dmd.convert_to_cgs()
tdyn = ne.evaluate("sqrt(3*pi/(32*G*(dens + dmd)))")
tdyn = tdyn*yt.units.s
return tdyn
def _sfr(field, data):
dens = data['density']
tdyn = data['total_dynamical_time']
Mu = np.float64(data.ds.parameters['Mu'])
dens.convert_to_cgs()
tdyn.convert_to_cgs()
sfr = ne.evaluate("where(dens*Mu/mh > 50, 0.01*dens/tdyn, 0)")
sfr = sfr*(yt.units.g/yt.units.cm**3/yt.units.s)
return sfr
ds.add_field('total_dynamical_time', _total_dynamical_time, units='s')
ds.add_field('sfr_density', _sfr, units='g/cm**3/s')
In [30]:
prj = yt.ProjectionPlot(ds, 2, 'sfr_density', width=(10, 'kpc'))
prj.set_unit('sfr_density', 'Msun/kpc**2/yr')
prj.set_zlim('sfr_density', 1e-1, 3e1)
Out[30]:
PhasePlot and ProfilePlot¶
In [31]:
from yt import PhasePlot
PhasePlot(ds.all_data(), 'density', 'temperature', 'cell_mass',
weight_field=None, fractional=True)
Out[31]:
In [33]:
max_val, max_loc = ds.find_max(('gas', 'density'))
print (max_val, max_loc)
In [34]:
PhasePlot(ds.sphere(max_loc, (3, 'kpc')), 'density', 'temperature', 'cell_mass',
weight_field=None, fractional=True)
Out[34]:
