yt.analysis_modules.halo_analysis.enzofof_merger_tree module

A very simple, purely-serial, merger tree script that knows how to parse FOF catalogs, either output by Enzo or output by yt’s FOF halo finder, and then compare parent/child relationships.

class yt.analysis_modules.halo_analysis.enzofof_merger_tree.EnzoFOFMergerBranch(tree, output_num, halo_id, max_children, min_relation=0.25)[source]

Bases: object

class yt.analysis_modules.halo_analysis.enzofof_merger_tree.EnzoFOFMergerTree(zrange=None, cycle_range=None, output=False, load_saved=False, save_filename='merger_tree.cpkl', external_FOF=True, FOF_directory='FOF')[source]

Bases: object

Calculates the parentage relationships for halos for a series of outputs, using the framework provided in enzofof_merger_tree.

Parameters:
  • zrange (tuple) – This is the redshift range (min, max) to calculate the merger tree. E.g. (0, 2) for z=2 to z=0
  • cycle_range (tuple, optional) – This is the cycle number range (min, max) to calculate the merger tree. If both zrange and cycle_number given, ignore zrange.
  • output (bool, optional) – If provided, both .cpkl and .txt files containing the parentage relationships will be output.
  • load_saved (bool, optional) – Flag to load previously saved parental relationships
  • save_filename (str, optional) – Filename to save parental relationships
  • external_FOF (bool, optional) – Are we building a tree from outputs generated by an external FOF program, or an FOF internal to yt?
  • FOF_directory (str, optional) – Directory where FOF files are located, note that the files must be named according to the syntax: groups_DDDDD.txt for internal yt outputs, and groups_DDDDD.dat for external FOF outputs. where DDDDD are digits representing the equivalent cycle number. e.g. groups_00000.txt

Examples

>>> mt = EnzoFOFMergerTree()    # by default it grabs every DD in FOF dir
>>> mt.build_tree(0)  # Create tree for halo 0
>>> mt.print_tree()
>>> mt.write_dot()
build_tree(halonum, min_particles=0, max_children=1e+20)[source]

Builds a merger tree, starting at the last output.

Parameters:
  • halonum (int) – Halo number in the last output to analyze.
  • min_particles (int, optional) – Minimum number of particles of halos in tree.
  • max_children (int, optional) – Maximum number of child halos each leaf can have.
clear_data()[source]

Deletes previous merger tree, but keeps parentage relationships.

filter_small_halos(lvl, min_particles)[source]
find_outputs(zrange, cycle_range, output)[source]
generate_tree(min_particles, max_children)[source]
get_massive_progenitors(halonum, min_relation=0.25)[source]

Returns a list of the most massive progenitor halos.

This routine walks down the tree, following the most massive progenitor on each node.

Parameters:halonum (int) – Halo number at the last output to trace.
Returns:output – Dictionary of redshifts, cycle numbers, and halo numbers of the most massive progenitor. keys = {redshift, cycle, halonum}
Return type:dict
load_tree(filename)[source]
print_tree()[source]

Prints the merger tree to stdout.

run_merger_tree(output)[source]
save_halo_evolution(filename)[source]

Saves as an HDF5 file the relevant details about a halo over the course of its evolution following the most massive progenitor to have given it the bulk of its particles. It stores info from the FOF_groups file: location, mass, id, etc.

save_tree(filename)[source]
select_cycles(cycle_range)[source]

Takes an existing tree and pares it to only include a subset of cycles. Useful in paring a loaded tree.

select_redshifts(zrange)[source]

Takes an existing tree and pares it to only include a subset of redshifts. Useful in paring a loaded tree.

write_dot(filename=None)[source]

Writes merger tree to a GraphViz or image file.

Parameters:filename (str, optional) – Filename to write the GraphViz file. Default will be tree_halo%05i.gv, which is a text file in the GraphViz format. If filename is an image (e.g. “MergerTree.png”) the output will be in the appropriate image format made by calling GraphViz automatically. See GraphViz (e.g. “dot -v”) for a list of available output formats.
class yt.analysis_modules.halo_analysis.enzofof_merger_tree.HaloCatalog(output_id, cache=True, external_FOF=True, FOF_directory='FOF')[source]

Bases: object

A catalog of halos, parsed from EnzoFOF outputs.

This class will read in catalogs output by the Enzo FOF halo finder and make available their positions, radii, etc. Enzo FOF was provided starting with 2.0, and can be run either inline (with the correct options) or as a postprocessing step using the -F command line option. This class is mostly useful when calculating a merger tree, and when the particle IDs for members of a given halo are output as well.

Parameters:
  • output_id (int) – This is the integer output id of the halo catalog to parse and load.
  • cache (bool) – Should we store, in between accesses, the particle IDs? If set to true, the correct particle files must exist.
  • external_FOF (bool, optional) – Are we building a tree from outputs generated by an external FOF program, or an FOF internal to yt?
  • FOF_directory (str, optional) – Directory where FOF files are located
cache = None
calculate_parentage_fractions(other_catalog, radius=0.1)[source]
parse_halo_catalog_external()[source]
parse_halo_catalog_internal()[source]

This parser works on the files output directly out of yt’s internal halo_finder. The parse_halo_catalog_external works with an external version of FOF.

Examples

>>> ds = load("DD0000/DD0000")
>>> halo_list = FOFHaloFinder(ds)
>>> halo_list.write_out("FOF/groups_00000.txt")
>>> halos_COM = parse_halo_catalog_internal()
read_particle_ids(halo_id)[source]
class yt.analysis_modules.halo_analysis.enzofof_merger_tree.HaloParticleList(halo_id, position, particle_ids)[source]

Bases: object

find_nearest(other_tree, radius=0.1)[source]
find_relative_parentage(child)[source]
class yt.analysis_modules.halo_analysis.enzofof_merger_tree.MaxLengthDict(*args, **kwargs)[source]

Bases: dict

clear() → None. Remove all items from D.
copy() → a shallow copy of D
fromkeys()

Returns a new dict with keys from iterable and values equal to value.

get(k[, d]) → D[k] if k in D, else d. d defaults to None.
items() → a set-like object providing a view on D's items
keys() → a set-like object providing a view on D's keys
pop(k[, d]) → v, remove specified key and return the corresponding value.

If key is not found, d is returned if given, otherwise KeyError is raised

popitem() → (k, v), remove and return some (key, value) pair as a

2-tuple; but raise KeyError if D is empty.

setdefault(k[, d]) → D.get(k,d), also set D[k]=d if k not in D
update([E, ]**F) → None. Update D from dict/iterable E and F.

If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]

values() → an object providing a view on D's values
yt.analysis_modules.halo_analysis.enzofof_merger_tree.find_halo_relationships(output1_id, output2_id, output_basename=None, radius=0.1, external_FOF=True, FOF_directory='FOF')[source]

Calculate the parentage and child relationships between two EnzoFOF halo catalogs.

This function performs a very simple merger tree calculation between two sets of halos. For every halo in the second halo catalog, it looks to the first halo catalog to find the parents by looking at particle IDs. The particle IDs from the child halos are identified in potential parents, and then both percent-of-parent and percent-to-child values are recorded.

Note that this works with catalogs constructed by Enzo’s FOF halo when used in external_FOF=True mode, whereas it will work with catalogs constructed by yt using external_FOF=False mode.

Parameters:
  • output1_id (int) – This is the integer output id of the (first) halo catalog to parse and load.
  • output2_id (int) – This is the integer output id of the (second) halo catalog to parse and load.
  • output_basename (string) – If provided, both .cpkl and .txt files containing the parentage relationships will be output.
  • radius (float, default to 0.10) – In absolute units, the radius to examine when guessing possible parent/child relationships. If this value is too small, you will miss possible relationships.
  • FOF_directory (str, optional) – Directory where FOF files are located
Returns:

pfrac – This is a dict of dicts. The first key is the parent halo id, the second is the child halo id. The values are the percent contributed from parent to child and the percent of a child that came from the parent.

Return type:

dict

yt.analysis_modules.halo_analysis.enzofof_merger_tree.grab_FOF_halo_info_internal(filename, halo_id)[source]

Finds a specific halo’s information in the FOF group output information and pass relevant parameters to caller.

yt.analysis_modules.halo_analysis.enzofof_merger_tree.plot_halo_evolution(filename, halo_id, x_quantity='cycle', y_quantity='mass', x_log=False, y_log=True, FOF_directory='FOF')[source]

Once you have generated a file using the EnzoFOFMergerTree.save_halo_evolution function, this is a simple way of plotting the evolution in the quantities of that halo over its lifetime.

Parameters:
  • filename (str) – The filename to which you saved the hdf5 data from save_halo_evolution
  • halo_id (int) – The halo in ‘filename’ that you want to follow
  • x_quantity (str, optional) –

    The quantity that you want to plot as the x_coord. Valid options are:

    • cycle
    • mass
    • fraction
    • halo_id
    • redshift
    • dense_x
    • dense_y
    • dense_z
    • COM_x
    • COM_y
    • COM_z
    • COM_vx
    • COM_vy
    • COM_vz
  • y_quantity (str, optional) – The quantity that you want to plot as the y_coord.
  • x_log (bool, optional) – Do you want the x-axis to be in log or linear?
  • y_log (bool, optional) – Do you want the y-axis to be in log or linear?
  • FOF_directory (str, optional) – Directory where FOF files (and hdf file) are located

Examples

>>> # generates mass history plots for the 20 most massive halos at t_fin.
>>> ts = DatasetSeries.from_filenames("DD????/DD????")
>>> # long step--must run FOF on each DD, but saves outputs for later use
>>> for ds in ts:
...     halo_list = FOFHaloFinder(ds)
...     i = int(ds.basename[2:])
...     halo_list.write_out("FOF/groups_%05i.txt" % i)
...     halo_list.write_particle_lists("FOF/particles_%05i" % i)
...
>>> mt = EnzoFOFMergerTree(external_FOF=False)
>>> for i in range(20):
...     mt.build_tree(i)
...     mt.save_halo_evolution('halos.h5')
...
>>> for i in range(20):
...     plot_halo_evolution('halos.h5', i)