[yt-svn] commit/yt: MatthewTurk: Merged in brittonsmith/yt/yt-3.0 (pull request #833)
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Sat Apr 19 04:59:13 PDT 2014
1 new commit in yt:
https://bitbucket.org/yt_analysis/yt/commits/6a60b95383d0/
Changeset: 6a60b95383d0
Branch: yt-3.0
User: MatthewTurk
Date: 2014-04-19 13:59:06
Summary: Merged in brittonsmith/yt/yt-3.0 (pull request #833)
Adding Halo Finder Registry and Subhalo Filter
Affected #: 6 files
diff -r 7de0dd178c22e32b843a58494d5bd6b1ffc2112c -r 6a60b95383d0ec478a532354de02dac6a62451d8 yt/analysis_modules/halo_analysis/api.py
--- a/yt/analysis_modules/halo_analysis/api.py
+++ b/yt/analysis_modules/halo_analysis/api.py
@@ -20,6 +20,9 @@
from .halo_callbacks import \
add_callback
+from .halo_finding_methods import \
+ add_finding_method
+
from .halo_filters import \
add_filter
diff -r 7de0dd178c22e32b843a58494d5bd6b1ffc2112c -r 6a60b95383d0ec478a532354de02dac6a62451d8 yt/analysis_modules/halo_analysis/finding_methods.py
--- a/yt/analysis_modules/halo_analysis/finding_methods.py
+++ /dev/null
@@ -1,20 +0,0 @@
-"""
-Halo Finding methods
-
-
-
-"""
-
-#-----------------------------------------------------------------------------
-# Copyright (c) 2013, yt Development Team.
-#
-# Distributed under the terms of the Modified BSD License.
-#
-# The full license is in the file COPYING.txt, distributed with this software.
-#-----------------------------------------------------------------------------
-
-from .operator_registry import \
- hf_registry
-
-class HaloFindingMethod(object):
- pass
diff -r 7de0dd178c22e32b843a58494d5bd6b1ffc2112c -r 6a60b95383d0ec478a532354de02dac6a62451d8 yt/analysis_modules/halo_analysis/halo_catalog.py
--- a/yt/analysis_modules/halo_analysis/halo_catalog.py
+++ b/yt/analysis_modules/halo_analysis/halo_catalog.py
@@ -30,20 +30,9 @@
from .operator_registry import \
callback_registry, \
filter_registry, \
- hf_registry, \
+ finding_method_registry, \
quantity_registry
-from yt.analysis_modules.halo_finding.halo_objects import \
- FOFHaloFinder, HOPHaloFinder
-from yt.frontends.halo_catalogs.halo_catalog.data_structures import \
- HaloCatalogDataset
-from yt.frontends.stream.data_structures import \
- load_particles
-from yt.frontends.halo_catalogs.rockstar.data_structures import \
- RockstarDataset
-from yt.analysis_modules.halo_finding.rockstar.api import \
- RockstarHaloFinder
-
class HaloCatalog(ParallelAnalysisInterface):
r"""Create a HaloCatalog: an object that allows for the creation and association
of data with a set of halo objects.
@@ -103,7 +92,7 @@
See Also
--------
- add_callback, add_filter, add_quantity
+ add_callback, add_filter, add_finding_method, add_quantity
"""
@@ -113,7 +102,6 @@
ParallelAnalysisInterface.__init__(self)
self.halos_pf = halos_pf
self.data_pf = data_pf
- self.finder_method = finder_method
self.output_dir = ensure_dir(output_dir)
if os.path.basename(self.output_dir) != ".":
self.output_prefix = os.path.basename(self.output_dir)
@@ -133,6 +121,10 @@
data_source = data_pf.h.all_data()
self.data_source = data_source
+ if finder_method is not None:
+ finder_method = finding_method_registry.find(finder_method)
+ self.finder_method = finder_method
+
# all of the analysis actions to be performed: callbacks, filters, and quantities
self.actions = []
# fields to be written to the halo catalog
@@ -358,16 +350,14 @@
if self.halos_pf is None:
# Find the halos and make a dataset of them
- particles_pf = self.find_halos()
+ self.halos_pf = self.finder_method(self.data_pf)
# Assign pf and data sources appropriately
- self.halos_pf = particles_pf
- self.data_source = particles_pf.all_data()
+ self.data_source = self.halos_pf.all_data()
# Add all of the default quantities that all halos must have
self.add_default_quantities('all')
-
my_index = np.argsort(self.data_source["particle_identifier"])
for i in parallel_objects(my_index, njobs=njobs, dynamic=dynamic):
new_halo = Halo(self)
@@ -400,80 +390,6 @@
if save_catalog:
self.save_catalog()
- def find_halos(self):
-
- finder_method = (self.finder_method).lower()
-
- if finder_method == "hop":
- halo_list = HOPHaloFinder(self.data_pf)
- halos_pf = self._parse_old_halo_list(halo_list)
-
- elif finder_method == "fof":
- halo_list = FOFHaloFinder(self.data_pf)
- halos_pf = self._parse_old_halo_list(halo_list)
-
- elif finder_method == 'rockstar':
- rh = RockstarHaloFinder(self.data_pf,
- outbase='{0}/rockstar_halos'.format(self.output_prefix))
- rh.run()
- halos_pf = RockstarDataset('{0}/rockstar_halos/halos_0.0.bin'.format(self.output_prefix))
- halos_pf.create_field_info()
- else:
- raise RuntimeError("finder_method must be 'fof', 'hop', or 'rockstar'")
-
- for attr in ["current_redshift", "current_time",
- "domain_dimensions",
- "cosmological_simulation", "omega_lambda",
- "omega_matter", "hubble_constant"]:
- attr_val = getattr(self.data_pf, attr)
- setattr(halos_pf, attr, attr_val)
- halos_pf.current_time = halos_pf.current_time.in_cgs()
-
- return halos_pf
-
- def _parse_old_halo_list(self, halo_list):
-
-
- data_pf = self.data_pf
- num_halos = len(halo_list)
-
- # Set up fields that we want to pull from identified halos and their units
- new_fields = ['particle_identifier', 'particle_mass', 'particle_position_x',
- 'particle_position_y','particle_position_z',
- 'virial_radius']
- new_units = [ '', 'g', 'cm', 'cm','cm','cm']
-
- # Set up a dictionary based on those fields
- # with empty arrays where we will fill in their values
- halo_properties = { f : (np.zeros(num_halos),unit) \
- for f, unit in zip(new_fields,new_units)}
-
- # Iterate through the halos pulling out their positions and virial quantities
- # and filling in the properties dictionary
- for i,halo in enumerate(halo_list):
- halo_properties['particle_identifier'][0][i] = i
- halo_properties['particle_mass'][0][i] = halo.virial_mass().in_cgs()
- halo_properties['virial_radius'][0][i] = halo.virial_radius().in_cgs()
-
- com = halo.center_of_mass().in_cgs()
- halo_properties['particle_position_x'][0][i] = com[0]
- halo_properties['particle_position_y'][0][i] = com[1]
- halo_properties['particle_position_z'][0][i] = com[2]
-
- # Define a bounding box based on original data pf
- bbox = np.array([data_pf.domain_left_edge.in_cgs(),
- data_pf.domain_right_edge.in_cgs()]).T
-
- # Create a pf with the halos as particles
- particle_pf = load_particles(halo_properties,
- bbox=bbox, length_unit = 1, mass_unit=1)
-
- # Create the field info dictionary so we can reference those fields
- particle_pf.create_field_info()
-
- return particle_pf
-
-
def save_catalog(self):
"Write out hdf5 file with all halo quantities."
@@ -513,4 +429,3 @@
self.add_quantity("particle_position_z", field_type=field_type)
self.add_quantity("virial_radius", field_type=field_type)
-
diff -r 7de0dd178c22e32b843a58494d5bd6b1ffc2112c -r 6a60b95383d0ec478a532354de02dac6a62451d8 yt/analysis_modules/halo_analysis/halo_filters.py
--- a/yt/analysis_modules/halo_analysis/halo_filters.py
+++ b/yt/analysis_modules/halo_analysis/halo_filters.py
@@ -13,6 +13,10 @@
# The full license is in the file COPYING.txt, distributed with this software.
#-----------------------------------------------------------------------------
+import numpy as np
+
+from yt.utilities.spatial import KDTree
+
from .halo_callbacks import HaloCallback
from .operator_registry import filter_registry
@@ -58,3 +62,44 @@
return eval("%s %s %s" % (h_value, operator, value))
add_filter("quantity_value", quantity_value)
+
+def _not_subhalo(halo, field_type="halos"):
+ """
+ Only return true if this halo is not a subhalo.
+
+ This is used for halo finders such as Rockstar that output parent
+ and subhalos together.
+ """
+
+ if not hasattr(halo.halo_catalog, "parent_dict"):
+ halo.halo_catalog.parent_dict = \
+ create_parent_dict(halo.halo_catalog.data_source, ptype=field_type)
+ return halo.halo_catalog.parent_dict[int(halo.quantities["particle_identifier"])] == -1
+add_filter("not_subhalo", _not_subhalo)
+
+def create_parent_dict(data_source, ptype="halos"):
+ """
+ Create a dictionary of halo parents to allow for filtering of subhalos.
+
+ For a pair of halos whose distance is smaller than the radius of at least
+ one of the halos, the parent is defined as the halo with the larger radius.
+ Parent halos (halos with no parents of their own) have parent index values of -1.
+ """
+ pos = np.rollaxis(
+ np.array([data_source[ptype, "particle_position_x"].in_units("Mpc"),
+ data_source[ptype, "particle_position_y"].in_units("Mpc"),
+ data_source[ptype, "particle_position_z"].in_units("Mpc")]), 1)
+ rad = data_source[ptype, "virial_radius"].in_units("Mpc").to_ndarray()
+ ids = data_source[ptype, "particle_identifier"].to_ndarray().astype("int")
+ parents = -1 * np.ones_like(ids, dtype="int")
+ my_tree = KDTree(pos)
+
+ for i in xrange(ids.size):
+ neighbors = np.array(
+ my_tree.query_ball_point(pos[i], rad[i], p=2))
+ if neighbors.size > 1:
+ parents[neighbors] = ids[neighbors[np.argmax(rad[neighbors])]]
+
+ parents[ids == parents] = -1
+ parent_dict = dict(zip(ids, parents))
+ return parent_dict
diff -r 7de0dd178c22e32b843a58494d5bd6b1ffc2112c -r 6a60b95383d0ec478a532354de02dac6a62451d8 yt/analysis_modules/halo_analysis/halo_finding_methods.py
--- /dev/null
+++ b/yt/analysis_modules/halo_analysis/halo_finding_methods.py
@@ -0,0 +1,132 @@
+"""
+Halo Finding methods
+
+
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+import numpy as np
+
+from yt.analysis_modules.halo_finding.halo_objects import \
+ FOFHaloFinder, HOPHaloFinder
+from yt.frontends.halo_catalogs.halo_catalog.data_structures import \
+ HaloCatalogDataset
+from yt.frontends.stream.data_structures import \
+ load_particles
+
+from .operator_registry import \
+ finding_method_registry
+
+def add_finding_method(name, function):
+ finding_method_registry[name] = HaloFindingMethod(function)
+
+class HaloFindingMethod(object):
+ r"""
+ A halo finding method is a callback that performs halo finding on a
+ dataset and returns a new dataset that is the loaded halo finder output.
+ """
+ def __init__(self, function, args=None, kwargs=None):
+ self.function = function
+ self.args = args
+ if self.args is None: self.args = []
+ self.kwargs = kwargs
+ if self.kwargs is None: self.kwargs = {}
+
+ def __call__(self, ds):
+ return self.function(ds, *self.args, **self.kwargs)
+
+def _hop_method(pf):
+ r"""
+ Run the Hop halo finding method.
+ """
+
+ halo_list = HOPHaloFinder(pf)
+ halos_pf = _parse_old_halo_list(pf, halo_list)
+ return halos_pf
+add_finding_method("hop", _hop_method)
+
+def _fof_method(pf):
+ r"""
+ Run the FoF halo finding method.
+ """
+
+ halo_list = FOFHaloFinder(pf)
+ halos_pf = _parse_old_halo_list(pf, halo_list)
+ return halos_pf
+add_finding_method("fof", _fof_method)
+
+def _rockstar_method(pf):
+ r"""
+ Run the Rockstar halo finding method.
+ """
+
+ from yt.frontends.halo_catalogs.rockstar.data_structures import \
+ RockstarDataset
+ from yt.analysis_modules.halo_finding.rockstar.api import \
+ RockstarHaloFinder
+
+ rh = RockstarHaloFinder(pf)
+ rh.run()
+ halos_pf = RockstarDataset("rockstar_halos/halos_0.0.bin")
+ halos_pf.create_field_info()
+ return halos_pf
+add_finding_method("rockstar", _rockstar_method)
+
+def _parse_old_halo_list(data_pf, halo_list):
+ r"""
+ Convert the halo list into a loaded dataset.
+ """
+
+ num_halos = len(halo_list)
+
+ # Set up fields that we want to pull from identified halos and their units
+ new_fields = ['particle_identifier', 'particle_mass', 'particle_position_x',
+ 'particle_position_y','particle_position_z',
+ 'virial_radius']
+ new_units = [ '', 'g', 'cm', 'cm','cm','cm']
+
+ # Set up a dictionary based on those fields
+ # with empty arrays where we will fill in their values
+ halo_properties = { f : (np.zeros(num_halos),unit) \
+ for f, unit in zip(new_fields,new_units)}
+
+ # Iterate through the halos pulling out their positions and virial quantities
+ # and filling in the properties dictionary
+ for i,halo in enumerate(halo_list):
+ halo_properties['particle_identifier'][0][i] = i
+ halo_properties['particle_mass'][0][i] = halo.virial_mass().in_cgs()
+ halo_properties['virial_radius'][0][i] = halo.virial_radius().in_cgs()
+
+ com = halo.center_of_mass().in_cgs()
+ halo_properties['particle_position_x'][0][i] = com[0]
+ halo_properties['particle_position_y'][0][i] = com[1]
+ halo_properties['particle_position_z'][0][i] = com[2]
+
+ # Define a bounding box based on original data pf
+ bbox = np.array([data_pf.domain_left_edge.in_cgs(),
+ data_pf.domain_right_edge.in_cgs()]).T
+
+ # Create a pf with the halos as particles
+ particle_pf = load_particles(halo_properties,
+ bbox=bbox, length_unit = 1, mass_unit=1)
+
+ # Create the field info dictionary so we can reference those fields
+ particle_pf.create_field_info()
+
+ for attr in ["current_redshift", "current_time",
+ "domain_dimensions",
+ "cosmological_simulation", "omega_lambda",
+ "omega_matter", "hubble_constant"]:
+ attr_val = getattr(data_pf, attr)
+ setattr(particle_pf, attr, attr_val)
+ particle_pf.current_time = particle_pf.current_time.in_cgs()
+
+ return particle_pf
diff -r 7de0dd178c22e32b843a58494d5bd6b1ffc2112c -r 6a60b95383d0ec478a532354de02dac6a62451d8 yt/analysis_modules/halo_analysis/operator_registry.py
--- a/yt/analysis_modules/halo_analysis/operator_registry.py
+++ b/yt/analysis_modules/halo_analysis/operator_registry.py
@@ -27,5 +27,5 @@
callback_registry = OperatorRegistry()
filter_registry = OperatorRegistry()
-hf_registry = OperatorRegistry()
+finding_method_registry = OperatorRegistry()
quantity_registry = OperatorRegistry()
Repository URL: https://bitbucket.org/yt_analysis/yt/
--
This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.
More information about the yt-svn
mailing list