Creating Derived Fields

One of the more powerful means of extending yt is through the usage of derived fields. These are fields that describe a value at each cell in a simulation.

Defining a New Field

So once a new field has been conceived of, the best way to create it is to construct a function that performs an array operation – operating on a collection of data, neutral to its size, shape, and type. (All fields should be provided as 64-bit floats.)

A simple example of this is the pressure field, which demonstrates the ease of this approach.

def _Pressure(field, data):
    return (data.pf["Gamma"] - 1.0) * \
           data["Density"] * data["ThermalEnergy"]

Note that we do a couple different things here. We access the “Gamma” parameter from the parameter file, we access the “Density” field and we access the “ThermalEnergy” field. “ThermalEnergy” is, in fact, another derived field! (“ThermalEnergy” deals with the distinction in storage of energy between dual energy formalism and non-DEF.) We don’t do any loops, we don’t do any type-checking, we can simply multiply the three items together.

Once we’ve defined our function, we need to notify yt that the field is available. The add_field() function is the means of doing this; it has a number of fairly specific parameters that can be passed in, but here we’ll only look at the most basic ones needed for a simple scalar baryon field.

add_field("Pressure", function=_Pressure, units=r"\rm{dyne}/\rm{cm}^{2}")

We feed it the name of the field, the name of the function, and the units. Note that the units parameter is a “raw” string, with some LaTeX-style formatting – Matplotlib actually has a MathText rendering engine, so if you include LaTeX it will be rendered appropriately.

We suggest that you name the function that creates a derived field with the intended field name prefixed by a single underscore, as in the _Pressure example above.

Note one last thing about this definition; we do not do unit conversion. All of the fields fed into the field are pre-supposed to be in CGS. If the field does not need any constants applied after that, you are done. If it does, you should define a second function that applies the proper multiple in order to return the desired units and use the argument convert_function to add_field to point to it.

If you find yourself using the same custom-defined fields over and over, you should put them in your plugins file as described in The Plugin File.

Some More Complicated Examples

But what if we want to do some more fancy stuff? Here’s an example of getting parameters from the data object and using those to define the field; specifically, here we obtain the center and height_vector parameters and use those to define an angle of declination of a point with respect to a disk.

def _DiskAngle(field, data):
    # We make both r_vec and h_vec into unit vectors
    center = data.get_field_parameter("center")
    r_vec = na.array([data["x"] - center[0],
                      data["y"] - center[1],
                      data["z"] - center[2]])
    r_vec = r_vec/na.sqrt((r_vec**2.0).sum(axis=0))
    h_vec = na.array(data.get_field_parameter("height_vector"))
    dp = r_vec[0,:] * h_vec[0] \
       + r_vec[1,:] * h_vec[1] \
       + r_vec[2,:] * h_vec[2]
    return na.arccos(dp)
add_field("DiskAngle", take_log=False,
          validators=[ValidateParameter("height_vector"),
                      ValidateParameter("center")],
          display_field=False)

Note that we have added a few parameters below the main function; we specify that we do not wish to display this field as logged, that we require both height_vector and center to be present in a given data object we wish to calculate this for, and we say that it should not be displayed in a drop-down box of fields to display. This is done through the parameter validators, which accepts a list of FieldValidator objects. These objects define the way in which the field is generated, and when it is able to be created. In this case, we mandate that parameters center and height_vector are set before creating the field. These are set via set_field_parameter(), which can be called on any object that has fields.

We can also define vector fields.

def _SpecificAngularMomentum(field, data):
    if data.has_field_parameter("bulk_velocity"):
        bv = data.get_field_parameter("bulk_velocity")
    else: bv = na.zeros(3, dtype='float64')
    xv = data["x-velocity"] - bv[0]
    yv = data["y-velocity"] - bv[1]
    zv = data["z-velocity"] - bv[2]
    center = data.get_field_parameter('center')
    coords = na.array([data['x'],data['y'],data['z']], dtype='float64')
    new_shape = tuple([3] + [1]*(len(coords.shape)-1))
    r_vec = coords - na.reshape(center,new_shape)
    v_vec = na.array([xv,yv,zv], dtype='float64')
    return na.cross(r_vec, v_vec, axis=0)
def _convertSpecificAngularMomentum(data):
    return data.convert("cm")
add_field("SpecificAngularMomentum",
          convert_function=_convertSpecificAngularMomentum, vector_field=True,
          units=r"\rm{cm}^2/\rm{s}", validators=[ValidateParameter('center')])

Here we define the SpecificAngularMomentum field, optionally taking a bulk_velocity, and returning a vector field that needs conversion by the function _convertSpecificAngularMomentum.

Field Options

The arguments to add_field() are passed on to the constructor of DerivedField. add_field() takes care of finding the arguments function and convert_function if it can, however. There are a number of options available, but the only mandatory ones are name and possibly function.

name
This is the name of the field – how you refer to it. For instance, Pressure or H2I_Fraction.
function
This is a function handle that defines the field
convert_function
This is the function that converts the field to CGS. All inputs to this function are mandated to already be in CGS.
units
This is a mathtext (LaTeX-like) string that describes the units.
projected_units
This is a mathtext (LaTeX-like) string that describes the units if the field has been projected without a weighting.
display_name
This is a name used in the plots
take_log
This is True or False and describes whether the field should be logged when plotted.
particle_type
Is this field a particle field?
validators
(Advanced) This is a list of FieldValidator objects, for instance to mandate spatial data.
vector_field
(Advanced) Is this field more than one value per cell?
display_field
(Advanced) Should this field appear in the dropdown box in Reason?
not_in_all
(Advanced) If this is True, the field may not be in all the grids.
projection_conversion
(Advanced) Which unit should we multiply by in a projection?

How Do Units Work?

Everything is done under the assumption that all of the native Enzo fields that yt knows about are converted to cgs before being handed to any processing routines.

Which Enzo Fields Does yt Know About?

  • Density
  • Temperature
  • Gas Energy
  • Total Energy
  • [xyz]-velocity
  • Species fields: HI, HII, Electron, HeI, HeII, HeIII, HM, H2I, H2II, DI, DII, HDI
  • Particle mass, velocity,

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