Using and Manipulating Objects and Fields

To generate standard plots, objects rarely need to be directly constructed. However, for detailed data inspection as well as hand-crafted derived data, objects can be exceptionally useful and even necessary.

Accessing Fields in Objects

yt utilizes load-on-demand objects to represent physical regions in space. (See Object Methodology.) Data objects in yt all respect the following protocol for accessing data:

my_object["Density"]

where "Density" can be any field name. The full list of objects is available in Available Objects, and information about how to create an object can be found in Creating 3D Datatypes. The field is returned as a single, flattened array without spatial information. The best mechanism for manipulating spatial data is the CoveringGridBase object.

The full list of fields that are available can be found as a property of the Hierarchy or Static Output object that you wish to access. This property is calculated every time the object is instantiated. The full list of fields that have been identified in the output file, which need no processing (besides unit conversion) are in the property field_list and the full list of potentially-accessible derived fields (see Derived Fields and Derived Quantities) is available in the property derived_field_list. You can see these by examining the two properties:

pf = load("my_data")
print pf.h.field_list
print pf.h.derived_field_list

When a field is added, it is added to a container that hangs off of the parameter file, as well. All of the field creation options (Field Options) are accessible through this object:

pf = load("my_data")
print pf.h.field_info["Pressure"].units

This is a fast way to examine the units of a given field, and additionally you can use yt.lagos.DerivedField.get_source() to get the source code:

field = pf.h.field_info["Pressure"]
print field.get_source()

Available Objects

Objects are instantiated by direct access of a hierarchy. Each of the objects that can be generated by a hierarchy are in fact fully-fledged data objects respecting the standard protocol for interaction.

The following objects are available, all of which hang off of the hierarchy object. To access them, you would do something like this (as for a region):

from yt.mods import *
pf = load("RedshiftOutput0005")
reg = pf.h.region([0.5, 0.5, 0.5], [0.0, 0.0, 0.0], [1.0, 1.0, 1.0])
covering_grid(self, level, left_edge, right_edge, dims, fields=None, pf=None, num_ghost_zones=0, use_pbar=True, **field_parameters):
(This is a proxy for yt.lagos.AMRCoveringGridBase.) The data object returned will consider grids up to level in generating fixed resolution data between left_edge and right_edge that is dims (3-values) on a side.
cutting(self, normal, center, fields=None, node_name=None, **field_parameters):
(This is a proxy for yt.lagos.AMRCuttingPlaneBase.) The Cutting Plane slices at an oblique angle, where we use the normal vector and the center to define the viewing plane. The ‘up’ direction is guessed at automatically.
disk(self, center, normal, radius, height, fields=None, pf=None, **field_parameters):
(This is a proxy for yt.lagos.AMRCylinderBase.) By providing a center, a normal, a radius and a height we can define a cylinder of any proportion. Only cells whose centers are within the cylinder will be selected.
extracted_region(self, base_region, indices, force_refresh=True, **field_parameters):
(This is a proxy for yt.lagos.ExtractedRegionBase.) Returns an instance of AMR3DData, or prepares one. Usually only used as a base class. Note that center is supplied, but only used for fields and quantities that require it.
grid(self, id, filename=None, hierarchy=None):
(This is a proxy for yt.lagos.EnzoGridBase.) Returns an instance of EnzoGrid with id, associated with filename and hierarchy.
grid_collection(self, center, grid_list, fields=None, pf=None, **field_parameters):
(This is a proxy for yt.lagos.AMRGridCollectionBase.) By selecting an arbitrary grid_list, we can act on those grids. Child cells are not returned.
ortho_ray(self, axis, coords, fields=None, pf=None, **field_parameters):
(This is a proxy for yt.lagos.AMROrthoRayBase.) Dimensionality is reduced to one, and an ordered list of points at an (x,y) tuple along axis are available.
periodic_region(self, center, left_edge, right_edge, fields=None, pf=None, **field_parameters):
(This is a proxy for yt.lagos.AMRPeriodicRegionBase.) We create an object with a set of three left_edge coordinates, three right_edge coordinates, and a center that need not be the center.
periodic_region_strict(self, center, left_edge, right_edge, fields=None, pf=None, **field_parameters):
(This is a proxy for yt.lagos.AMRPeriodicRegionStrictBase.) We create an object with a set of three left_edge coordinates, three right_edge coordinates, and a center that need not be the center.
proj(self, axis, field, weight_field=None, max_level=None, center=None, pf=None, source=None, node_name=None, field_cuts=None, serialize=True, **field_parameters):
(This is a proxy for yt.lagos.AMRProjBase.) AMRProj is a projection of a field along an axis. The field can have an associated weight_field, in which case the values are multiplied by a weight before being summed, and then divided by the sum of that weight.
ray(self, start_point, end_point, fields=None, pf=None, **field_parameters):
(This is a proxy for yt.lagos.AMRRayBase.) We accept a start point and an end point and then get all the data between those two.
region(self, center, left_edge, right_edge, fields=None, pf=None, **field_parameters):
(This is a proxy for yt.lagos.AMRRegionBase.) We create an object with a set of three left_edge coordinates, three right_edge coordinates, and a center that need not be the center.
region_strict(self, center, left_edge, right_edge, fields=None, pf=None, **field_parameters):
(This is a proxy for yt.lagos.AMRRegionStrictBase.) We create an object with a set of three left_edge coordinates, three right_edge coordinates, and a center that need not be the center.
slice(self, axis, coord, fields=None, center=None, pf=None, node_name=False, **field_parameters):
(This is a proxy for yt.lagos.AMRSliceBase.) Slice along axisHow do I specify an axis?, at the coordinate coord. Optionally supply fields.
smoothed_covering_grid(self, *args, **field_parameters):
(This is a proxy for yt.lagos.AMRSmoothedCoveringGridBase.) The data object returned will consider grids up to level in generating fixed resolution data between left_edge and right_edge that is dims (3-values) on a side.
sphere(self, center, radius, fields=None, pf=None, **field_parameters):
(This is a proxy for yt.lagos.AMRSphereBase.) The most famous of all the data objects, we define it via a center and a radius.

Storing and Loading Objects

Often, when operating interactively or via the scripting interface, it is convenient to save an object or multiple objects out to disk and then restart the calculation later. Personally, I found this most useful when dealing with identification of clumps and contours (see Cookbook for a recipe on how to find clumps and the API documentation for both ContourFinder and Clump) where the identification step can be quite time-consuming, but the analysis may be relatively fast.

Typically, the save and load operations are used on 3D data objects. yt has a separate set of serialization operations for 2D objects such as projections.

yt will save out 3D objects to disk under the presupposition that the construction of the objects is the difficult part, rather than the generation of the data – this means that you can save out an object as a description of how to recreate it in space, but not the actual data arrays affiliated with that object. The information that is saved includes the parameter file off of which the object “hangs.” It is this piece of information that is the most difficult; the object, when reloaded, must be able to reconstruct a parameter file from whatever limited information it has in the save file.

To do this, yt is able to identify parameter files based on a “hash” generated from the base file name, the “CurrentTimeIdentifier”, and the simulation time. These three characteristics should never be changed outside of a simulation, they are independent of the file location on disk, and in conjunction they should be uniquely identifying. (This process is all done in fido via ParameterFileStore.)

To save an object, you can either save it in the .yt file affiliated with the hierarchy (What are all these .yt files?) or as a standalone file. For instance, using save_object() we can save a sphere.

from yt.mods import *
pf = load("my_data")
sp = pf.h.sphere([0.5, 0.5, 0.5], 10.0/pf['kpc'])

pf.h.save_object(sp, "sphere_to_analyze_later")

In a later session, we can load it using save_object():

from yt.mods import *

pf = load("my_data")
sphere_to_analyze = pf.h.load_object("sphere_to_analyze_later")

Additionally, if we want to store the object independent of the .yt file, we can save the object directly:

from yt.mods import *

pf = load("my_data")
sp = pf.h.sphere([0.5, 0.5, 0.5], 10.0/pf['kpc'])

sp.save_object("my_sphere", "my_storage_file.cpkl")

This will store the object as my_sphere in the file my_storage_file.cpkl, which will be created or accessed using the standard python module shelve. Note that if a filename is not supplied, it will be saved via the hierarchy, as above.

To re-load an object saved this way, you can use the shelve module directly:

from yt.mods import *
import shelve

obj_file = shelve.open("my_storage_file.cpkl")
pf, obj = obj_file["my_sphere"]

Note here that this behaves slightly differently than above – we do not need to load the parameter file ourselves, as the load process actually does that for us! Additionally, we can store multiple objects in a single shelve file, so we have to call the sphere by name.

Note

It’s also possible to use the standard cPickle module for loading and storing objects – so in theory you could even save a list of objects!

This method works for clumps, as well, and the entire clump hierarchy will be stored and restored upon load.

Comments

Feel free to leave comments! If you've got a GMail account, you can use https://www.google.com/accounts/o8/id as your OpenID URL.
comments powered by Disqus

Table Of Contents

Previous topic

Command Line Tool

Next topic

Examining and Manipulating Particles

This Page