diff --git a/docs/install.rst b/docs/install.rst
index a6b937276..f9e1c8cd3 100644
--- a/docs/install.rst
+++ b/docs/install.rst
@@ -14,10 +14,14 @@ The simplest way to install the latest release of the code is with pip. Before i
pip install halotools
-Alternatively, you can install using conda::
+Alternatively, you can install using the astropy channel on conda::
conda install -c astropy halotools
+or via conda-forge:
+
+ conda install -c conda-forge halotools
+
Either pip or conda will install the latest official release of the code.
If instead you want the latest master branch,
you will need to build the code from source following the instructions in the next section.
@@ -58,13 +62,13 @@ The first step is to clone the halotools repository::
Installing one of the official releases
------------------------------------------
-All official releases of the code are tagged with their version name, e.g., v0.5.
+All official releases of the code are tagged with their version name, e.g., v0.7.
To install a particular release::
- git checkout v0.5
+ git checkout v0.7
python setup.py install
-This will install the v0.5 release of the code. Other official release versions (e.g., v0.1) can be installed similarly.
+This will install the v0.7 release of the code. Other official release versions (e.g., v0.5) can be installed similarly.
Installing the most recent master branch
------------------------------------------
@@ -113,7 +117,7 @@ The h5py package is used for fast I/O of large simulated datasets.
If you did not use pip, then you should be aware of the following strict requirements:
-- `Python `_: 2.7.x or 3.x
+- `Python `_: 3.7.x
- `Numpy `_: 1.9 or later
@@ -121,7 +125,7 @@ If you did not use pip, then you should be aware of the following strict require
- `Cython `_: 0.23 or later
-- `Astropy`_: 1.0 or later
+- `Astropy`_: 4.0 or later
- `BeautifulSoup `_: For crawling the web for halo catalogs.
diff --git a/docs/installing_halotools_with_virtualenv.rst b/docs/installing_halotools_with_virtualenv.rst
index ebb9baab8..bb4c31dc0 100644
--- a/docs/installing_halotools_with_virtualenv.rst
+++ b/docs/installing_halotools_with_virtualenv.rst
@@ -12,11 +12,11 @@ create a virtual environment that will automatically have compatible versions of
By installing into a virtual environment, you will not change any of the
packages that are already installed system-wide on your machine. In the example below, we will use conda to create a virtual environment with all the dependencies handled automatically::
- conda create -n halotools_env python=3.7 astropy numpy scipy h5py requests beautifulsoup4 cython
+ conda create -n ht07 python=3.7 halotools=0.7 h5py ipython jupyter matplotlib
In order to activate this environment::
- source activate halotools_env
+ source activate ht07
Then install halotools into this environment::
@@ -24,7 +24,7 @@ Then install halotools into this environment::
Or, alternatively, you can install the latest master branch by following the :ref:`install_from_source` instructions.
-Any additional packages you install into the ``halotools_env`` virtual environment will not impact your system-wide environment. Whenever you want to do science involving Halotools,
+Any additional packages you install into the ``ht07`` virtual environment will not impact your system-wide environment. Whenever you want to do science involving Halotools,
just activate the environment and import the code. When you are done
and wish to return to your normal system environment::
diff --git a/halotools/empirical_models/composite_models/hod_models/tests/test_zheng07.py b/halotools/empirical_models/composite_models/hod_models/tests/test_zheng07.py
index 20f495f05..2cddcf469 100644
--- a/halotools/empirical_models/composite_models/hod_models/tests/test_zheng07.py
+++ b/halotools/empirical_models/composite_models/hod_models/tests/test_zheng07.py
@@ -10,18 +10,22 @@
from .....custom_exceptions import HalotoolsError
from .....sim_manager import FakeSim
-__all__ = ('test_zheng07_composite1', 'test_zheng07_composite2')
+__all__ = ("test_zheng07_composite1", "test_zheng07_composite2")
def test_zheng07_composite1():
""" Ensure that the ``zheng07`` pre-built model does not accept non-Zheng07Cens.
"""
- model1 = PrebuiltHodModelFactory('zheng07')
- model2 = PrebuiltHodModelFactory('zheng07', modulate_with_cenocc=True)
+ model1 = PrebuiltHodModelFactory("zheng07")
+ model2 = PrebuiltHodModelFactory("zheng07", modulate_with_cenocc=True)
with pytest.raises(HalotoolsError) as err:
- __ = PrebuiltHodModelFactory('zheng07', modulate_with_cenocc=True, cenocc_model=0)
- substr = "Do not pass in the ``cenocc_model`` keyword to ``zheng07_model_dictionary``"
+ __ = PrebuiltHodModelFactory(
+ "zheng07", modulate_with_cenocc=True, cenocc_model=0
+ )
+ substr = (
+ "Do not pass in the ``cenocc_model`` keyword to ``zheng07_model_dictionary``"
+ )
assert substr in err.value.args[0]
@@ -31,31 +35,39 @@ def test_zheng07_composite2():
"""
class MyCenModel(OccupationComponent):
-
def __init__(self, threshold):
- OccupationComponent.__init__(self, gal_type='centrals',
- threshold=threshold, upper_occupation_bound=1.)
-
- self.param_dict['new_cen_param'] = 0.5
+ OccupationComponent.__init__(
+ self,
+ gal_type="centrals",
+ threshold=threshold,
+ upper_occupation_bound=1.0,
+ )
+
+ self.param_dict["new_cen_param"] = 0.5
self._suppress_repeated_param_warning = True
def mean_occupation(self, **kwargs):
- halo_table = kwargs['table']
- result = np.zeros(len(halo_table)) + self.param_dict['new_cen_param']
+ halo_table = kwargs["table"]
+ result = np.zeros(len(halo_table)) + self.param_dict["new_cen_param"]
return result
centrals_occupation = MyCenModel(threshold=-20)
- satellites_occupation = Zheng07Sats(threshold=-20,
- modulate_with_cenocc=True, cenocc_model=centrals_occupation,
- _suppress_repeated_param_warning=True)
+ satellites_occupation = Zheng07Sats(
+ threshold=-20,
+ modulate_with_cenocc=True,
+ cenocc_model=centrals_occupation,
+ _suppress_repeated_param_warning=True,
+ )
centrals_profile = TrivialPhaseSpace()
satellites_profile = NFWPhaseSpace()
- model_dict = ({'centrals_occupation': centrals_occupation,
- 'centrals_profile': centrals_profile,
- 'satellites_occupation': satellites_occupation,
- 'satellites_profile': satellites_profile})
+ model_dict = {
+ "centrals_occupation": centrals_occupation,
+ "centrals_profile": centrals_profile,
+ "satellites_occupation": satellites_occupation,
+ "satellites_profile": satellites_profile,
+ }
composite_model = HodModelFactory(**model_dict)
@@ -67,27 +79,42 @@ def test_modulate_with_cenocc1():
""" Regression test for Issue #646. Verify that the ``modulate_with_cenocc``
keyword argument is actually passed to the satellite occupation component.
"""
- model1 = PrebuiltHodModelFactory('zheng07', modulate_with_cenocc=True)
- assert model1.model_dictionary[u'satellites_occupation'].modulate_with_cenocc is True
+ model1 = PrebuiltHodModelFactory("zheng07", modulate_with_cenocc=True)
+ assert (
+ model1.model_dictionary[u"satellites_occupation"].modulate_with_cenocc is True
+ )
def test_modulate_with_cenocc2():
""" Regression test for Issue #646. Verify that the ``modulate_with_cenocc``
keyword argument results in behavior that is properly modified by the centrals.
"""
- model = PrebuiltHodModelFactory('zheng07', modulate_with_cenocc=True)
+ model = PrebuiltHodModelFactory("zheng07", modulate_with_cenocc=True)
test_mass = 1e12
ncen1 = model.mean_occupation_satellites(prim_haloprop=test_mass)
- model.param_dict['logMmin'] *= 1.5
+ model.param_dict["logMmin"] *= 1.5
ncen2 = model.mean_occupation_satellites(prim_haloprop=test_mass)
assert ncen1 != ncen2
def test_host_centric_distance():
- model = PrebuiltHodModelFactory('zheng07')
- halocat = FakeSim(seed=43)
- model.populate_mock(halocat, seed=43)
+ model = PrebuiltHodModelFactory("zheng07")
+ halocat = FakeSim(seed=43)
+ model.populate_mock(halocat, seed=43)
- msg = "host_centric_distance key was never mapped during mock population"
- assert np.any(model.mock.galaxy_table['host_centric_distance'] > 0), msg
+ msg = "host_centric_distance key was never mapped during mock population"
+ assert np.any(model.mock.galaxy_table["host_centric_distance"] > 0), msg
+
+def test_mass_definition_flexibility():
+ """ Regression test for Issue #993.
+ """
+ model = PrebuiltHodModelFactory(
+ "zheng07", mdef="200m", halo_boundary_key="halo_radius_arbitrary"
+ )
+ halocat = FakeSim(seed=43)
+ halocat.halo_table["halo_m200m"] = np.copy(halocat.halo_table["halo_mvir"])
+ halocat.halo_table["halo_radius_arbitrary"] = np.copy(
+ halocat.halo_table["halo_mvir"]
+ )
+ model.populate_mock(halocat, seed=43)
diff --git a/halotools/empirical_models/factories/hod_mock_factory.py b/halotools/empirical_models/factories/hod_mock_factory.py
index dbb31c421..1cab463ea 100644
--- a/halotools/empirical_models/factories/hod_mock_factory.py
+++ b/halotools/empirical_models/factories/hod_mock_factory.py
@@ -22,16 +22,20 @@
from ...custom_exceptions import HalotoolsError
-__all__ = ['HodMockFactory']
-__author__ = ['Andrew Hearin']
+__all__ = ["HodMockFactory"]
+__author__ = ["Andrew Hearin"]
-unavailable_haloprop_msg = ("Your model requires that the ``%s`` key appear in the halo catalog,\n"
- "but this column is not available in the catalog you attempted to populate.\n")
+unavailable_haloprop_msg = (
+ "Your model requires that the ``%s`` key appear in the halo catalog,\n"
+ "but this column is not available in the catalog you attempted to populate.\n"
+)
-missing_halo_upid_msg = ("All HOD-style models populate host halos with mock galaxies.\n"
+missing_halo_upid_msg = (
+ "All HOD-style models populate host halos with mock galaxies.\n"
"The way Halotools distinguishes host halos from subhalos is via the ``halo_upid`` column,\n"
"with halo_upid = -1 for host halos and !=-1 for subhalos.\n"
- "The halo catalog you passed to the HodMockFactory does not have the ``halo_upid`` column.\n")
+ "The halo catalog you passed to the HodMockFactory does not have the ``halo_upid`` column.\n"
+)
class HodMockFactory(MockFactory):
@@ -51,8 +55,12 @@ class HodMockFactory(MockFactory):
"""
- def __init__(self, Num_ptcl_requirement=sim_defaults.Num_ptcl_requirement,
- halo_mass_column_key='halo_mvir', **kwargs):
+ def __init__(
+ self,
+ Num_ptcl_requirement=sim_defaults.Num_ptcl_requirement,
+ halo_mass_column_key="halo_mvir",
+ **kwargs
+ ):
"""
Parameters
----------
@@ -86,7 +94,7 @@ def __init__(self, Num_ptcl_requirement=sim_defaults.Num_ptcl_requirement,
MockFactory.__init__(self, **kwargs)
- halocat = kwargs['halocat']
+ halocat = kwargs["halocat"]
self.Num_ptcl_requirement = Num_ptcl_requirement
self.halo_mass_column_key = halo_mass_column_key
@@ -115,22 +123,24 @@ def preprocess_halo_catalog(self, halocat):
"""
try:
- assert 'halo_upid' in list(halocat.halo_table.keys())
+ assert "halo_upid" in list(halocat.halo_table.keys())
except AssertionError:
raise HalotoolsError(missing_halo_upid_msg)
# Make cuts on halo catalog #
# Select host halos only, since this is an HOD-style model
halo_table, subhalo_table = SampleSelector.host_halo_selection(
- table=halocat.halo_table, return_subhalos=True)
+ table=halocat.halo_table, return_subhalos=True
+ )
# make a (possibly trivial) completeness cut
- cutoff_mvir = self.Num_ptcl_requirement*self.particle_mass
+ cutoff_mvir = self.Num_ptcl_requirement * self.particle_mass
mass_cut = halo_table[self.halo_mass_column_key] > cutoff_mvir
max_column_value = np.max(halo_table[self.halo_mass_column_key])
halo_table = halo_table[mass_cut]
if len(halo_table) == 0:
- msg = ("During the pre-processing phase of HOD mock-making \n"
+ msg = (
+ "During the pre-processing phase of HOD mock-making \n"
"controlled by the `preprocess_halo_catalog` method of the HodMockFactory,\n"
"a cut on halo completeness is made. The column used in this cut\n"
"is determined by the value of the halo_mass_column_key string \n"
@@ -142,15 +152,23 @@ def preprocess_halo_catalog(self, halocat):
"the processing of the halo catalog,\n"
"for example, an incorrect column number and/or dtype.\n"
"Alternatively, the value of Num_ptcl_requirement = %.2e \n"
- "that was passed to the HodMockFactory constructor could be the problem.\n")
- raise HalotoolsError(msg % (
- self.halo_mass_column_key, max_column_value, cutoff_mvir,
- self.particle_mass, self.Num_ptcl_requirement))
+ "that was passed to the HodMockFactory constructor could be the problem.\n"
+ )
+ raise HalotoolsError(
+ msg
+ % (
+ self.halo_mass_column_key,
+ max_column_value,
+ cutoff_mvir,
+ self.particle_mass,
+ self.Num_ptcl_requirement,
+ )
+ )
# If any component model needs the subhalo_table, bind it to the mock object
for component_model in self.model.model_dictionary.values():
try:
- f = getattr(component_model, 'preprocess_subhalo_table')
+ f = getattr(component_model, "preprocess_subhalo_table")
halo_table, self.subhalo_table = f(halo_table, subhalo_table)
except AttributeError:
pass
@@ -289,17 +307,17 @@ def populate(self, seed=None, **kwargs):
# The _testing_mode keyword is for unit-testing only
# it has been intentionally left out of the docstring
try:
- self._testing_mode = kwargs['_testing_mode']
+ self._testing_mode = kwargs["_testing_mode"]
except KeyError:
self._testing_mode = False
try:
- self.enforce_PBC = kwargs['enforce_PBC']
+ self.enforce_PBC = kwargs["enforce_PBC"]
except KeyError:
self.enforce_PBC = True
try:
- masking_function = kwargs['masking_function']
+ masking_function = kwargs["masking_function"]
mask = masking_function(self._orig_halo_table)
self.halo_table = self._orig_halo_table[mask]
except:
@@ -318,21 +336,23 @@ def populate(self, seed=None, **kwargs):
# For the gal_type_slice indices of
# the pre-allocated array self.gal_type,
# set each string-type entry equal to the gal_type string
- self.galaxy_table['gal_type'][gal_type_slice] = (
- np.repeat(gal_type, self._total_abundance[gal_type], axis=0))
+ self.galaxy_table["gal_type"][gal_type_slice] = np.repeat(
+ gal_type, self._total_abundance[gal_type], axis=0
+ )
# Store all other relevant host halo properties into their
# appropriate pre-allocated array
for halocatkey in self.additional_haloprops:
self.galaxy_table[halocatkey][gal_type_slice] = np.repeat(
- self.halo_table[halocatkey], self._occupation[gal_type], axis=0)
+ self.halo_table[halocatkey], self._occupation[gal_type], axis=0
+ )
- self.galaxy_table['x'] = self.galaxy_table['halo_x']
- self.galaxy_table['y'] = self.galaxy_table['halo_y']
- self.galaxy_table['z'] = self.galaxy_table['halo_z']
- self.galaxy_table['vx'] = self.galaxy_table['halo_vx']
- self.galaxy_table['vy'] = self.galaxy_table['halo_vy']
- self.galaxy_table['vz'] = self.galaxy_table['halo_vz']
+ self.galaxy_table["x"] = self.galaxy_table["halo_x"]
+ self.galaxy_table["y"] = self.galaxy_table["halo_y"]
+ self.galaxy_table["z"] = self.galaxy_table["halo_z"]
+ self.galaxy_table["vx"] = self.galaxy_table["halo_vx"]
+ self.galaxy_table["vy"] = self.galaxy_table["halo_vy"]
+ self.galaxy_table["vz"] = self.galaxy_table["halo_vz"]
for method in self._remaining_methods_to_call:
func = getattr(self.model, method)
@@ -346,28 +366,37 @@ def populate(self, seed=None, **kwargs):
func(table=self.galaxy_table[gal_type_slice], seed=seed, **d)
if self.enforce_PBC is True:
- self.galaxy_table['x'], self.galaxy_table['vx'] = (
- model_helpers.enforce_periodicity_of_box(
- self.galaxy_table['x'], self.Lbox[0],
- velocity=self.galaxy_table['vx'],
- check_multiple_box_lengths=self._testing_mode)
- )
-
- self.galaxy_table['y'], self.galaxy_table['vy'] = (
- model_helpers.enforce_periodicity_of_box(
- self.galaxy_table['y'], self.Lbox[1],
- velocity=self.galaxy_table['vy'],
- check_multiple_box_lengths=self._testing_mode)
- )
-
- self.galaxy_table['z'], self.galaxy_table['vz'] = (
- model_helpers.enforce_periodicity_of_box(
- self.galaxy_table['z'], self.Lbox[2],
- velocity=self.galaxy_table['vz'],
- check_multiple_box_lengths=self._testing_mode)
- )
-
- if hasattr(self.model, 'galaxy_selection_func'):
+ (
+ self.galaxy_table["x"],
+ self.galaxy_table["vx"],
+ ) = model_helpers.enforce_periodicity_of_box(
+ self.galaxy_table["x"],
+ self.Lbox[0],
+ velocity=self.galaxy_table["vx"],
+ check_multiple_box_lengths=self._testing_mode,
+ )
+
+ (
+ self.galaxy_table["y"],
+ self.galaxy_table["vy"],
+ ) = model_helpers.enforce_periodicity_of_box(
+ self.galaxy_table["y"],
+ self.Lbox[1],
+ velocity=self.galaxy_table["vy"],
+ check_multiple_box_lengths=self._testing_mode,
+ )
+
+ (
+ self.galaxy_table["z"],
+ self.galaxy_table["vz"],
+ ) = model_helpers.enforce_periodicity_of_box(
+ self.galaxy_table["z"],
+ self.Lbox[2],
+ velocity=self.galaxy_table["vz"],
+ check_multiple_box_lengths=self._testing_mode,
+ )
+
+ if hasattr(self.model, "galaxy_selection_func"):
mask = self.model.galaxy_selection_func(self.galaxy_table)
self.galaxy_table = self.galaxy_table[mask]
@@ -386,14 +415,16 @@ def allocate_memory(self, seed=None):
# We will keep track of the calling sequence with a list called _remaining_methods_to_call
# Each time a function in this list is called, we will remove that function from the list
# Mock generation will be complete when _remaining_methods_to_call is exhausted
- self._remaining_methods_to_call = copy(self.model._mock_generation_calling_sequence)
+ self._remaining_methods_to_call = copy(
+ self.model._mock_generation_calling_sequence
+ )
# Call all composite model methods that should be called prior to mc_occupation
# All such function calls must be applied to the table, since we do not yet know
# how much memory we need for the mock galaxy_table
galprops_assigned_to_halo_table = []
for func_name in self.model._mock_generation_calling_sequence:
- if 'mc_occupation' in func_name:
+ if "mc_occupation" in func_name:
# exit when we encounter a ``mc_occupation_`` function
break
else:
@@ -405,8 +436,12 @@ def allocate_memory(self, seed=None):
if seed is not None:
seed += 1
func(table=self.halo_table, seed=seed, **d)
- galprops_assigned_to_halo_table_by_func = func._galprop_dtypes_to_allocate.names
- galprops_assigned_to_halo_table.extend(galprops_assigned_to_halo_table_by_func)
+ galprops_assigned_to_halo_table_by_func = (
+ func._galprop_dtypes_to_allocate.names
+ )
+ galprops_assigned_to_halo_table.extend(
+ galprops_assigned_to_halo_table_by_func
+ )
self._remaining_methods_to_call.remove(func_name)
# Now update the list of additional_haloprops, if applicable
# This is necessary because each of the above function calls created new
@@ -414,8 +449,7 @@ def allocate_memory(self, seed=None):
# np.repeat inside mock.populate() so that mock galaxies inherit these newly-created columns
# Since there is already a loop over additional_haloprops inside mock.populate() that does this,
# then all we need to do is append to this list
- galprops_assigned_to_halo_table = list(set(
- galprops_assigned_to_halo_table))
+ galprops_assigned_to_halo_table = list(set(galprops_assigned_to_halo_table))
self.additional_haloprops.extend(galprops_assigned_to_halo_table)
self.additional_haloprops = list(set(self.additional_haloprops))
@@ -424,49 +458,51 @@ def allocate_memory(self, seed=None):
self._gal_type_indices = {}
for gal_type in self.gal_types:
- self.halo_table['halo_num_'+gal_type] = 0
+ self.halo_table["halo_num_" + gal_type] = 0
first_galaxy_index = 0
for gal_type in self.gal_types:
- occupation_func_name = 'mc_occupation_'+gal_type
+ occupation_func_name = "mc_occupation_" + gal_type
occupation_func = getattr(self.model, occupation_func_name)
# Call the component model to get a Monte Carlo
# realization of the abundance of gal_type galaxies
if seed is not None:
seed += 1
- self._occupation[gal_type] = occupation_func(table=self.halo_table, seed=seed)
- self.halo_table['halo_num_'+gal_type][:] = self._occupation[gal_type]
+ self._occupation[gal_type] = occupation_func(
+ table=self.halo_table, seed=seed
+ )
+ self.halo_table["halo_num_" + gal_type][:] = self._occupation[gal_type]
# Now use the above result to set up the indexing scheme
- self._total_abundance[gal_type] = (
- self._occupation[gal_type].sum()
- )
+ self._total_abundance[gal_type] = self._occupation[gal_type].sum()
last_galaxy_index = first_galaxy_index + self._total_abundance[gal_type]
# Build a bookkeeping device to keep track of
# which array elements pertain to which gal_type.
self._gal_type_indices[gal_type] = slice(
- first_galaxy_index, last_galaxy_index)
+ first_galaxy_index, last_galaxy_index
+ )
first_galaxy_index = last_galaxy_index
# Remove the mc_occupation function from the list of methods to call
self._remaining_methods_to_call.remove(occupation_func_name)
- self.additional_haloprops.append('halo_num_'+gal_type)
+ self.additional_haloprops.append("halo_num_" + gal_type)
self.Ngals = np.sum(list(self._total_abundance.values()))
# Allocate memory for all additional halo properties,
# including profile parameters of the halos such as 'conc_NFWmodel'
for halocatkey in self.additional_haloprops:
- self.galaxy_table[halocatkey] = np.zeros(self.Ngals,
- dtype=self.halo_table[halocatkey].dtype)
+ self.galaxy_table[halocatkey] = np.zeros(
+ self.Ngals, dtype=self.halo_table[halocatkey].dtype
+ )
# Separately allocate memory for the galaxy profile parameters
for galcatkey in self.model.halo_prof_param_keys:
- self.galaxy_table[galcatkey] = 0.
+ self.galaxy_table[galcatkey] = 0.0
for galcatkey in self.model.gal_prof_param_keys:
- self.galaxy_table[galcatkey] = 0.
+ self.galaxy_table[galcatkey] = 0.0
- self.galaxy_table['gal_type'] = np.zeros(self.Ngals, dtype=object)
+ self.galaxy_table["gal_type"] = np.zeros(self.Ngals, dtype=object)
dt = self.model._galprop_dtypes_to_allocate
for key in dt.names:
@@ -487,7 +523,7 @@ def estimate_ngals(self, seed=None):
# All such function calls must be applied to the table.
halo_table = Table(np.copy(self.halo_table))
for func_name in self.model._mock_generation_calling_sequence:
- if 'mc_occupation' in func_name:
+ if "mc_occupation" in func_name:
# exit when we encounter a ``mc_occupation_`` function
break
else:
@@ -498,7 +534,7 @@ def estimate_ngals(self, seed=None):
# realization of the abundance of galaxies for all gal_type.
ngals = 0
for gal_type in self.gal_types:
- occupation_func_name = 'mc_occupation_'+gal_type
+ occupation_func_name = "mc_occupation_" + gal_type
occupation_func = getattr(self.model, occupation_func_name)
if seed is not None:
seed += 1
diff --git a/halotools/empirical_models/phase_space_models/analytic_models/centrals/trivial_phase_space.py b/halotools/empirical_models/phase_space_models/analytic_models/centrals/trivial_phase_space.py
index 5cf0a899c..f5db0b494 100644
--- a/halotools/empirical_models/phase_space_models/analytic_models/centrals/trivial_phase_space.py
+++ b/halotools/empirical_models/phase_space_models/analytic_models/centrals/trivial_phase_space.py
@@ -11,8 +11,8 @@
from .....sim_manager import sim_defaults
-__author__ = ['Andrew Hearin']
-__all__ = ['TrivialPhaseSpace']
+__author__ = ["Andrew Hearin"]
+__all__ = ["TrivialPhaseSpace"]
class TrivialPhaseSpace(object):
@@ -22,11 +22,14 @@ class TrivialPhaseSpace(object):
:math:`P(\vec{x}_{\rm cen}, \vec{v}_{\rm cen}) = \delta^{(6)}(\vec{x}_{\rm halo}, \vec{v}_{\rm halo})`.
"""
- def __init__(self,
- cosmology=sim_defaults.default_cosmology,
- redshift=sim_defaults.default_redshift,
- mdef=model_defaults.halo_mass_definition,
- **kwargs):
+ def __init__(
+ self,
+ cosmology=sim_defaults.default_cosmology,
+ redshift=sim_defaults.default_redshift,
+ mdef=model_defaults.halo_mass_definition,
+ halo_boundary_key=None,
+ **kwargs
+ ):
"""
Parameters
----------
@@ -40,22 +43,31 @@ def __init__(self,
String specifying the halo mass definition, e.g., 'vir' or '200m'.
Default is set in `~halotools.empirical_models.model_defaults`.
"""
- self._mock_generation_calling_sequence = ['assign_phase_space']
- self._galprop_dtypes_to_allocate = np.dtype([
- ('x', 'f8'), ('y', 'f8'), ('z', 'f8'),
- ('vx', 'f8'), ('vy', 'f8'), ('vz', 'f8'),
- ])
+ self._mock_generation_calling_sequence = ["assign_phase_space"]
+ self._galprop_dtypes_to_allocate = np.dtype(
+ [
+ ("x", "f8"),
+ ("y", "f8"),
+ ("z", "f8"),
+ ("vx", "f8"),
+ ("vy", "f8"),
+ ("vz", "f8"),
+ ]
+ )
self.param_dict = {}
self.cosmology = cosmology
self.redshift = redshift
self.mdef = mdef
- self.halo_boundary_key = model_defaults.get_halo_boundary_key(self.mdef)
+ if halo_boundary_key is None:
+ self.halo_boundary_key = model_defaults.get_halo_boundary_key(self.mdef)
+ else:
+ self.halo_boundary_key = halo_boundary_key
def assign_phase_space(self, table, **kwargs):
r"""
"""
- phase_space_keys = ['x', 'y', 'z', 'vx', 'vy', 'vz']
+ phase_space_keys = ["x", "y", "z", "vx", "vy", "vz"]
for key in phase_space_keys:
- table[key][:] = table['halo_'+key][:]
+ table[key][:] = table["halo_" + key][:]
diff --git a/halotools/empirical_models/phase_space_models/analytic_models/profile_model_template.py b/halotools/empirical_models/phase_space_models/analytic_models/profile_model_template.py
index 54f35a988..16831938d 100644
--- a/halotools/empirical_models/phase_space_models/analytic_models/profile_model_template.py
+++ b/halotools/empirical_models/phase_space_models/analytic_models/profile_model_template.py
@@ -17,11 +17,11 @@
from ... import model_defaults
-newtonG = G.to(u.km*u.km*u.Mpc/(u.Msun*u.s*u.s))
+newtonG = G.to(u.km * u.km * u.Mpc / (u.Msun * u.s * u.s))
-__author__ = ['Andrew Hearin', 'Benedikt Diemer']
+__author__ = ["Andrew Hearin", "Benedikt Diemer"]
-__all__ = ['AnalyticDensityProf']
+__all__ = ["AnalyticDensityProf"]
@six.add_metaclass(ABCMeta)
@@ -42,7 +42,7 @@ class AnalyticDensityProf(object):
the profile is the necessary and sufficient ingredient.
"""
- def __init__(self, cosmology, redshift, mdef, **kwargs):
+ def __init__(self, cosmology, redshift, mdef, halo_boundary_key=None, **kwargs):
r"""
Parameters
-----------
@@ -55,6 +55,9 @@ def __init__(self, cosmology, redshift, mdef, **kwargs):
mdef: str
String specifying the halo mass definition, e.g., 'vir' or '200m'.
+ halo_boundary_key : str, optional
+ Default behavior is to use the column associated with the input mdef.
+
"""
self.cosmology = cosmology
self.redshift = redshift
@@ -63,9 +66,12 @@ def __init__(self, cosmology, redshift, mdef, **kwargs):
# The following four attributes are derived quantities from the above,
# so that self-consistency between them is ensured
self.density_threshold = halo_boundary_functions.density_threshold(
- cosmology=self.cosmology,
- redshift=self.redshift, mdef=self.mdef)
- self.halo_boundary_key = model_defaults.get_halo_boundary_key(self.mdef)
+ cosmology=self.cosmology, redshift=self.redshift, mdef=self.mdef
+ )
+ if halo_boundary_key is None:
+ self.halo_boundary_key = model_defaults.get_halo_boundary_key(self.mdef)
+ else:
+ self.halo_boundary_key = halo_boundary_key
self.prim_haloprop_key = model_defaults.get_halo_mass_key(self.mdef)
self.gal_prof_param_keys = []
@@ -157,11 +163,13 @@ def mass_density(self, radius, mass, *prof_params):
"""
halo_radius = self.halo_mass_to_halo_radius(mass)
- scaled_radius = radius/halo_radius
+ scaled_radius = radius / halo_radius
- dimensionless_mass = self.dimensionless_mass_density(scaled_radius, *prof_params)
+ dimensionless_mass = self.dimensionless_mass_density(
+ scaled_radius, *prof_params
+ )
- density = self.density_threshold*dimensionless_mass
+ density = self.density_threshold * dimensionless_mass
return density
def _enclosed_dimensionless_mass_integrand(self, scaled_radius, *prof_params):
@@ -183,8 +191,10 @@ def _enclosed_dimensionless_mass_integrand(self, scaled_radius, *prof_params):
integrand: array_like
function to be integrated to yield the amount of enclosed mass.
"""
- dimensionless_density = self.dimensionless_mass_density(scaled_radius, *prof_params)
- return dimensionless_density*4*np.pi*scaled_radius**2
+ dimensionless_density = self.dimensionless_mass_density(
+ scaled_radius, *prof_params
+ )
+ return dimensionless_density * 4 * np.pi * scaled_radius ** 2
def cumulative_mass_PDF(self, scaled_radius, *prof_params):
r"""
@@ -219,12 +229,20 @@ def cumulative_mass_PDF(self, scaled_radius, *prof_params):
for i in range(len(x)):
enclosed_mass[i], _ = quad_integration(
- self._enclosed_dimensionless_mass_integrand, 0., x[i], epsrel=1e-5,
- args=prof_params)
+ self._enclosed_dimensionless_mass_integrand,
+ 0.0,
+ x[i],
+ epsrel=1e-5,
+ args=prof_params,
+ )
total, _ = quad_integration(
- self._enclosed_dimensionless_mass_integrand, 0., 1.0, epsrel=1e-5,
- args=prof_params)
+ self._enclosed_dimensionless_mass_integrand,
+ 0.0,
+ 1.0,
+ epsrel=1e-5,
+ args=prof_params,
+ )
return enclosed_mass / total
@@ -259,7 +277,7 @@ def enclosed_mass(self, radius, total_mass, *prof_params):
"""
radius = np.atleast_1d(radius).astype(np.float64)
scaled_radius = radius / self.halo_mass_to_halo_radius(total_mass)
- mass = self.cumulative_mass_PDF(scaled_radius, *prof_params)*total_mass
+ mass = self.cumulative_mass_PDF(scaled_radius, *prof_params) * total_mass
return mass
@@ -289,7 +307,9 @@ def dimensionless_circular_velocity(self, scaled_radius, *prof_params):
See :ref:`halo_profile_definitions` for derivations and implementation details.
"""
- return np.sqrt(self.cumulative_mass_PDF(scaled_radius, *prof_params)/scaled_radius)
+ return np.sqrt(
+ self.cumulative_mass_PDF(scaled_radius, *prof_params) / scaled_radius
+ )
def virial_velocity(self, total_mass):
r""" The circular velocity evaluated at the halo boundary,
@@ -310,8 +330,9 @@ def virial_velocity(self, total_mass):
See :ref:`halo_profile_definitions` for derivations and implementation details.
"""
- return halo_boundary_functions.halo_mass_to_virial_velocity(total_mass,
- self.cosmology, self.redshift, self.mdef)
+ return halo_boundary_functions.halo_mass_to_virial_velocity(
+ total_mass, self.cosmology, self.redshift, self.mdef
+ )
def circular_velocity(self, radius, total_mass, *prof_params):
r"""
@@ -344,13 +365,14 @@ def circular_velocity(self, radius, total_mass, *prof_params):
halo_radius = self.halo_mass_to_halo_radius(total_mass)
scaled_radius = np.atleast_1d(radius) / halo_radius
return self.dimensionless_circular_velocity(
- scaled_radius, *prof_params)*self.virial_velocity(total_mass)
+ scaled_radius, *prof_params
+ ) * self.virial_velocity(total_mass)
def _vmax_helper(self, scaled_radius, *prof_params):
""" Helper function used to calculate `vmax` and `rmax`.
"""
encl = self.cumulative_mass_PDF(scaled_radius, *prof_params)
- return -1.*encl/scaled_radius
+ return -1.0 * encl / scaled_radius
def rmax(self, total_mass, *prof_params):
r""" Radius at which the halo attains its maximum circular velocity.
@@ -380,7 +402,7 @@ def rmax(self, total_mass, *prof_params):
result = scipy_minimize(self._vmax_helper, guess, args=prof_params)
- return result.x[0]*halo_radius
+ return result.x[0] * halo_radius
def vmax(self, total_mass, *prof_params):
r""" Maximum circular velocity of the halo profile.
@@ -408,7 +430,9 @@ def vmax(self, total_mass, *prof_params):
result = scipy_minimize(self._vmax_helper, guess, args=prof_params)
halo_radius = self.halo_mass_to_halo_radius(total_mass)
- return self.circular_velocity(result.x[0]*halo_radius, total_mass, *prof_params)
+ return self.circular_velocity(
+ result.x[0] * halo_radius, total_mass, *prof_params
+ )
def halo_mass_to_halo_radius(self, total_mass):
r"""
@@ -434,8 +458,8 @@ def halo_mass_to_halo_radius(self, total_mass):
"""
return halo_boundary_functions.halo_mass_to_halo_radius(
- total_mass, cosmology=self.cosmology,
- redshift=self.redshift, mdef=self.mdef)
+ total_mass, cosmology=self.cosmology, redshift=self.redshift, mdef=self.mdef
+ )
def halo_radius_to_halo_mass(self, radius):
r"""
@@ -461,5 +485,5 @@ def halo_radius_to_halo_mass(self, radius):
"""
return halo_boundary_functions.halo_radius_to_halo_mass(
- radius, cosmology=self.cosmology,
- redshift=self.redshift, mdef=self.mdef)
+ radius, cosmology=self.cosmology, redshift=self.redshift, mdef=self.mdef
+ )
diff --git a/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/biased_nfw_phase_space.py b/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/biased_nfw_phase_space.py
index fa7315cbb..6426b0b8b 100644
--- a/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/biased_nfw_phase_space.py
+++ b/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/biased_nfw_phase_space.py
@@ -14,16 +14,18 @@
from ..... import model_defaults
-__author__ = ('Andrew Hearin', )
-__all__ = ('BiasedNFWPhaseSpace', )
+__author__ = ("Andrew Hearin",)
+__all__ = ("BiasedNFWPhaseSpace",)
-lookup_table_performance_warning = ("You have selected {0} bins to digitize host halo concentration \n"
+lookup_table_performance_warning = (
+ "You have selected {0} bins to digitize host halo concentration \n"
"and {1} bins to digitize the galaxy bias parameter.\n"
"To populate mocks, the BiasedNFWPhaseSpace class builds a lookup table with shape ({0}, {1}, {2}),\n"
"one entry for every numerical solution to the Jeans equation.\n"
"Using this fine of a binning requires a long pre-computation of {3} integrals\n."
- "Make sure you actually need to use so many bins")
+ "Make sure you actually need to use so many bins"
+)
class BiasedNFWPhaseSpace(NFWPhaseSpace):
@@ -37,6 +39,7 @@ class BiasedNFWPhaseSpace(NFWPhaseSpace):
including descriptions of how the relevant equations are
implemented in the Halotools code base, see :ref:`nfw_profile_tutorial`.
"""
+
def __init__(self, profile_integration_tol=1e-5, **kwargs):
r"""
Parameters
@@ -68,6 +71,9 @@ def __init__(self, profile_integration_tol=1e-5, **kwargs):
String specifying the halo mass definition, e.g., 'vir' or '200m'.
Default is set in `~halotools.empirical_models.model_defaults`.
+ halo_boundary_key : str, optional
+ Default behavior is to use the column associated with the input mdef.
+
concentration_key : string, optional
Column name of the halo catalog storing NFW concentration.
@@ -142,20 +148,30 @@ def __init__(self, profile_integration_tol=1e-5, **kwargs):
NFWPhaseSpace.__init__(self, **kwargs)
self._profile_integration_tol = profile_integration_tol
- self.gal_prof_param_keys = ['conc_NFWmodel', 'conc_gal_bias']
+ self.gal_prof_param_keys = ["conc_NFWmodel", "conc_gal_bias"]
prof_lookup_bins = self._retrieve_prof_lookup_info(**kwargs)
self.setup_prof_lookup_tables(*prof_lookup_bins)
- self._mock_generation_calling_sequence = ['calculate_conc_gal_bias', 'assign_phase_space']
- self._methods_to_inherit = ['calculate_conc_gal_bias']
-
- self._galprop_dtypes_to_allocate = np.dtype([
- ('host_centric_distance', 'f8'),
- ('x', 'f8'), ('y', 'f8'), ('z', 'f8'),
- ('vx', 'f8'), ('vy', 'f8'), ('vz', 'f8'),
- ('conc_gal_bias', 'f8'), ('conc_galaxy', 'f8')
- ])
+ self._mock_generation_calling_sequence = [
+ "calculate_conc_gal_bias",
+ "assign_phase_space",
+ ]
+ self._methods_to_inherit = ["calculate_conc_gal_bias"]
+
+ self._galprop_dtypes_to_allocate = np.dtype(
+ [
+ ("host_centric_distance", "f8"),
+ ("x", "f8"),
+ ("y", "f8"),
+ ("z", "f8"),
+ ("vx", "f8"),
+ ("vy", "f8"),
+ ("vz", "f8"),
+ ("conc_gal_bias", "f8"),
+ ("conc_galaxy", "f8"),
+ ]
+ )
self.param_dict = self._initialize_conc_bias_param_dict(**kwargs)
@@ -163,37 +179,50 @@ def _initialize_conc_bias_param_dict(self, **kwargs):
r""" Set up the appropriate number of keys in the parameter dictionary
and give the keys standardized names.
"""
- if 'conc_gal_bias_logM_abscissa' in list(kwargs.keys()):
+ if "conc_gal_bias_logM_abscissa" in list(kwargs.keys()):
_conc_bias_logM_abscissa = np.atleast_1d(
- kwargs.get('conc_gal_bias_logM_abscissa')).astype('f4')
- d = ({'conc_gal_bias_param'+str(i): 1.
- for i in range(len(_conc_bias_logM_abscissa))})
- d2 = ({'conc_gal_bias_logM_abscissa_param'+str(i): float(logM)
- for i, logM in enumerate(_conc_bias_logM_abscissa)})
+ kwargs.get("conc_gal_bias_logM_abscissa")
+ ).astype("f4")
+ d = {
+ "conc_gal_bias_param" + str(i): 1.0
+ for i in range(len(_conc_bias_logM_abscissa))
+ }
+ d2 = {
+ "conc_gal_bias_logM_abscissa_param" + str(i): float(logM)
+ for i, logM in enumerate(_conc_bias_logM_abscissa)
+ }
self._num_conc_bias_params = len(_conc_bias_logM_abscissa)
for key, value in d2.items():
d[key] = value
return d
else:
- return {'conc_gal_bias': 1.}
+ return {"conc_gal_bias": 1.0}
def _retrieve_prof_lookup_info(self, **kwargs):
r""" Retrieve the arrays defining the lookup table control points
"""
- cmin, cmax = model_defaults.min_permitted_conc, model_defaults.max_permitted_conc
- dc = 1.
- npts_conc = int(np.round((cmax-cmin)/float(dc)))
+ cmin, cmax = (
+ model_defaults.min_permitted_conc,
+ model_defaults.max_permitted_conc,
+ )
+ dc = 1.0
+ npts_conc = int(np.round((cmax - cmin) / float(dc)))
default_conc_arr = np.linspace(cmin, cmax, npts_conc)
- conc_arr = kwargs.get('concentration_bins', default_conc_arr)
+ conc_arr = kwargs.get("concentration_bins", default_conc_arr)
- conc_gal_bias_arr = kwargs.get('conc_gal_bias_bins',
- model_defaults.default_conc_gal_bias_bins)
+ conc_gal_bias_arr = kwargs.get(
+ "conc_gal_bias_bins", model_defaults.default_conc_gal_bias_bins
+ )
npts_conc, npts_gal_bias = len(conc_arr), len(conc_gal_bias_arr)
- if npts_conc*npts_gal_bias > 300:
- args = (npts_conc, npts_gal_bias, model_defaults.Npts_radius_table,
- npts_conc*npts_gal_bias*model_defaults.Npts_radius_table)
+ if npts_conc * npts_gal_bias > 300:
+ args = (
+ npts_conc,
+ npts_gal_bias,
+ model_defaults.Npts_radius_table,
+ npts_conc * npts_gal_bias * model_defaults.Npts_radius_table,
+ )
warn(lookup_table_performance_warning.format(*args))
return [conc_arr, conc_gal_bias_arr]
@@ -235,7 +264,7 @@ def _clipped_galaxy_concentration(self, halo_conc, conc_gal_bias):
r""" Private method used to ensure that biased concentrations are still
bound by the min/max bounds permitted by the lookup tables.
"""
- gal_conc = conc_gal_bias*halo_conc
+ gal_conc = conc_gal_bias * halo_conc
try:
cmin = self._conc_NFWmodel_lookup_table_min
@@ -335,7 +364,9 @@ def cumulative_mass_PDF(self, scaled_radius, halo_conc):
"""
return NFWPhaseSpace.cumulative_mass_PDF(self, scaled_radius, halo_conc)
- def dimensionless_radial_velocity_dispersion(self, scaled_radius, halo_conc, conc_gal_bias):
+ def dimensionless_radial_velocity_dispersion(
+ self, scaled_radius, halo_conc, conc_gal_bias
+ ):
r"""
Analytical solution to the isotropic jeans equation for an NFW potential,
rendered dimensionless via scaling by the virial velocity.
@@ -361,8 +392,12 @@ def dimensionless_radial_velocity_dispersion(self, scaled_radius, halo_conc, con
The returned result has the same dimension as the input ``scaled_radius``.
"""
gal_conc = self._clipped_galaxy_concentration(halo_conc, conc_gal_bias)
- return biased_dimless_vrad_disp(scaled_radius, halo_conc, gal_conc,
- profile_integration_tol=self._profile_integration_tol)
+ return biased_dimless_vrad_disp(
+ scaled_radius,
+ halo_conc,
+ gal_conc,
+ profile_integration_tol=self._profile_integration_tol,
+ )
def radial_velocity_dispersion(self, radius, total_mass, halo_conc, conc_gal_bias):
r"""
@@ -394,13 +429,12 @@ def radial_velocity_dispersion(self, radius, total_mass, halo_conc, conc_gal_bia
"""
virial_velocities = self.virial_velocity(total_mass)
halo_radius = self.halo_mass_to_halo_radius(total_mass)
- scaled_radius = radius/halo_radius
+ scaled_radius = radius / halo_radius
- dimensionless_velocities = (
- self.dimensionless_radial_velocity_dispersion(
- scaled_radius, halo_conc, conc_gal_bias)
- )
- return dimensionless_velocities*virial_velocities
+ dimensionless_velocities = self.dimensionless_radial_velocity_dispersion(
+ scaled_radius, halo_conc, conc_gal_bias
+ )
+ return dimensionless_velocities * virial_velocities
def calculate_conc_gal_bias(self, seed=None, **kwargs):
r""" Calculate the ratio of the galaxy concentration to the halo concentration,
@@ -430,33 +464,39 @@ def calculate_conc_gal_bias(self, seed=None, **kwargs):
Ratio of the galaxy concentration to the halo concentration,
:math:`c_{\rm gal}/c_{\rm halo}`.
"""
- if 'table' in list(kwargs.keys()):
- table = kwargs['table']
+ if "table" in list(kwargs.keys()):
+ table = kwargs["table"]
mass = table[self.prim_haloprop_key]
- elif 'prim_haloprop' in list(kwargs.keys()):
- mass = np.atleast_1d(kwargs['prim_haloprop']).astype('f4')
+ elif "prim_haloprop" in list(kwargs.keys()):
+ mass = np.atleast_1d(kwargs["prim_haloprop"]).astype("f4")
else:
- msg = ("\nYou must pass either a ``table`` or ``prim_haloprop`` argument \n"
- "to the ``assign_conc_gal_bias`` function of the ``BiasedNFWPhaseSpace`` class.\n")
+ msg = (
+ "\nYou must pass either a ``table`` or ``prim_haloprop`` argument \n"
+ "to the ``assign_conc_gal_bias`` function of the ``BiasedNFWPhaseSpace`` class.\n"
+ )
raise KeyError(msg)
- if 'conc_gal_bias_logM_abscissa_param0' in self.param_dict.keys():
- abscissa_keys = list('conc_gal_bias_logM_abscissa_param'+str(i)
- for i in range(self._num_conc_bias_params))
+ if "conc_gal_bias_logM_abscissa_param0" in self.param_dict.keys():
+ abscissa_keys = list(
+ "conc_gal_bias_logM_abscissa_param" + str(i)
+ for i in range(self._num_conc_bias_params)
+ )
abscissa = [self.param_dict[key] for key in abscissa_keys]
- ordinates_keys = list('conc_gal_bias_param'+str(i)
- for i in range(self._num_conc_bias_params))
+ ordinates_keys = list(
+ "conc_gal_bias_param" + str(i)
+ for i in range(self._num_conc_bias_params)
+ )
ordinates = [self.param_dict[key] for key in ordinates_keys]
result = np.interp(np.log10(mass), abscissa, ordinates)
else:
- result = np.zeros_like(mass) + self.param_dict['conc_gal_bias']
+ result = np.zeros_like(mass) + self.param_dict["conc_gal_bias"]
- if 'table' in list(kwargs.keys()):
- table['conc_gal_bias'][:] = result
- halo_conc = table['conc_NFWmodel']
+ if "table" in list(kwargs.keys()):
+ table["conc_gal_bias"][:] = result
+ halo_conc = table["conc_NFWmodel"]
gal_conc = self._clipped_galaxy_concentration(halo_conc, result)
- table['conc_galaxy'][:] = gal_conc
+ table["conc_galaxy"][:] = gal_conc
else:
return result
diff --git a/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/nfw_phase_space.py b/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/nfw_phase_space.py
index a088ddda1..effb05b53 100644
--- a/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/nfw_phase_space.py
+++ b/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/nfw_phase_space.py
@@ -13,8 +13,8 @@
from ..... import model_defaults
-__author__ = ['Andrew Hearin']
-__all__ = ['NFWPhaseSpace']
+__author__ = ["Andrew Hearin"]
+__all__ = ["NFWPhaseSpace"]
class NFWPhaseSpace(NFWProfile, MonteCarloGalProf):
@@ -28,6 +28,7 @@ class NFWPhaseSpace(NFWProfile, MonteCarloGalProf):
including descriptions of how the relevant equations are
implemented in the Halotools code base, see :ref:`nfw_profile_tutorial`.
"""
+
def __init__(self, **kwargs):
r"""
Parameters
@@ -51,6 +52,9 @@ def __init__(self, **kwargs):
String specifying the halo mass definition, e.g., 'vir' or '200m'.
Default is set in `~halotools.empirical_models.model_defaults`.
+ halo_boundary_key : str, optional
+ Default behavior is to use the column associated with the input mdef.
+
concentration_key : string, optional
Column name of the halo catalog storing NFW concentration.
@@ -76,16 +80,19 @@ def __init__(self, **kwargs):
prof_lookup_args = self._retrieve_prof_lookup_info(**kwargs)
self.setup_prof_lookup_tables(*prof_lookup_args)
- self._mock_generation_calling_sequence = ['assign_phase_space']
+ self._mock_generation_calling_sequence = ["assign_phase_space"]
def _retrieve_prof_lookup_info(self, **kwargs):
r""" Retrieve the arrays defining the lookup table control points
"""
- cmin, cmax = model_defaults.min_permitted_conc, model_defaults.max_permitted_conc
- dc = 1.
- npts_conc = int(np.round((cmax-cmin)/float(dc)))
+ cmin, cmax = (
+ model_defaults.min_permitted_conc,
+ model_defaults.max_permitted_conc,
+ )
+ dc = 1.0
+ npts_conc = int(np.round((cmax - cmin) / float(dc)))
default_conc_arr = np.linspace(cmin, cmax, npts_conc)
- conc_arr = kwargs.get('concentration_bins', default_conc_arr)
+ conc_arr = kwargs.get("concentration_bins", default_conc_arr)
return [conc_arr]
def assign_phase_space(self, table, seed=None):
@@ -599,11 +606,12 @@ def radial_velocity_dispersion(self, radius, total_mass, halo_conc):
"""
virial_velocities = self.virial_velocity(total_mass)
halo_radius = self.halo_mass_to_halo_radius(total_mass)
- scaled_radius = radius/halo_radius
+ scaled_radius = radius / halo_radius
- dimensionless_velocities = (
- self.dimensionless_radial_velocity_dispersion(scaled_radius, halo_conc))
- return dimensionless_velocities*virial_velocities
+ dimensionless_velocities = self.dimensionless_radial_velocity_dispersion(
+ scaled_radius, halo_conc
+ )
+ return dimensionless_velocities * virial_velocities
def setup_prof_lookup_tables(self, *concentration_bins):
r"""
@@ -616,10 +624,12 @@ def setup_prof_lookup_tables(self, *concentration_bins):
"""
MonteCarloGalProf.setup_prof_lookup_tables(self, *concentration_bins)
- def build_lookup_tables(self,
- logrmin=model_defaults.default_lograd_min,
- logrmax=model_defaults.default_lograd_max,
- Npts_radius_table=model_defaults.Npts_radius_table):
+ def build_lookup_tables(
+ self,
+ logrmin=model_defaults.default_lograd_min,
+ logrmax=model_defaults.default_lograd_max,
+ Npts_radius_table=model_defaults.Npts_radius_table,
+ ):
r""" Method used to create a lookup table of the spatial and velocity radial profiles.
Parameters
@@ -715,7 +725,9 @@ def mc_halo_centric_pos(self, *concentration_array, **kwargs):
Length-Ngals array storing a Monte Carlo realization of the galaxy positions.
"""
- return MonteCarloGalProf.mc_halo_centric_pos(self, *concentration_array, **kwargs)
+ return MonteCarloGalProf.mc_halo_centric_pos(
+ self, *concentration_array, **kwargs
+ )
def mc_pos(self, *concentration_array, **kwargs):
r""" Method to generate random, three-dimensional positions of galaxies.
@@ -783,10 +795,13 @@ def _vrad_disp_from_lookup(self, scaled_radius, *concentration_array, **kwargs):
scaled by the size of the halo's virial velocity.
"""
- return MonteCarloGalProf._vrad_disp_from_lookup(self,
- scaled_radius, *concentration_array, **kwargs)
+ return MonteCarloGalProf._vrad_disp_from_lookup(
+ self, scaled_radius, *concentration_array, **kwargs
+ )
- def mc_radial_velocity(self, scaled_radius, total_mass, *concentration_array, **kwargs):
+ def mc_radial_velocity(
+ self, scaled_radius, total_mass, *concentration_array, **kwargs
+ ):
r"""
Method returns a Monte Carlo realization of radial velocities drawn from Gaussians
with a width determined by the solution to the isotropic Jeans equation.
@@ -813,8 +828,9 @@ def mc_radial_velocity(self, scaled_radius, total_mass, *concentration_array, **
Array of radial velocities drawn from Gaussians with a width determined by the
solution to the Jeans equation.
"""
- return MonteCarloGalProf.mc_radial_velocity(self,
- scaled_radius, total_mass, *concentration_array, **kwargs)
+ return MonteCarloGalProf.mc_radial_velocity(
+ self, scaled_radius, total_mass, *concentration_array, **kwargs
+ )
def mc_vel(self, table, seed=None):
r""" Method assigns a Monte Carlo realization of the Jeans velocity
@@ -834,8 +850,9 @@ def mc_vel(self, table, seed=None):
"""
MonteCarloGalProf.mc_vel(self, table, seed=seed)
- def mc_generate_nfw_phase_space_points(self, Ngals=int(1e4), conc=5, mass=1e12,
- verbose=True, seed=None):
+ def mc_generate_nfw_phase_space_points(
+ self, Ngals=int(1e4), conc=5, mass=1e12, verbose=True, seed=None
+ ):
r""" Return a Monte Carlo realization of points
in the phase space of an NFW halo in isotropic Jeans equilibrium.
@@ -913,10 +930,11 @@ def mc_generate_nfw_phase_space_points(self, Ngals=int(1e4), conc=5, mass=1e12,
rvir = NFWProfile.halo_mass_to_halo_radius(self, total_mass=m)
- x, y, z = MonteCarloGalProf.mc_halo_centric_pos(self, c,
- halo_radius=rvir, seed=seed)
- r = np.sqrt(x**2 + y**2 + z**2)
- scaled_radius = r/rvir
+ x, y, z = MonteCarloGalProf.mc_halo_centric_pos(
+ self, c, halo_radius=rvir, seed=seed
+ )
+ r = np.sqrt(x ** 2 + y ** 2 + z ** 2)
+ scaled_radius = r / rvir
if seed is not None:
seed += 1
@@ -932,16 +950,25 @@ def mc_generate_nfw_phase_space_points(self, Ngals=int(1e4), conc=5, mass=1e12,
yrel, vyrel = _relative_positions_and_velocities(y, 0, v1=vy, v2=0)
zrel, vzrel = _relative_positions_and_velocities(z, 0, v1=vz, v2=0)
- vrad = (xrel*vxrel + yrel*vyrel + zrel*vzrel)/r
-
- t = Table({'x': x, 'y': y, 'z': z,
- 'vx': vx, 'vy': vy, 'vz': vz,
- 'radial_position': r, 'radial_velocity': vrad})
+ vrad = (xrel * vxrel + yrel * vyrel + zrel * vzrel) / r
+
+ t = Table(
+ {
+ "x": x,
+ "y": y,
+ "z": z,
+ "vx": vx,
+ "vy": vy,
+ "vz": vz,
+ "radial_position": r,
+ "radial_velocity": vrad,
+ }
+ )
return t
-def _sign_pbc(x1, x2, period=None, equality_fill_val=0., return_pbc_correction=False):
+def _sign_pbc(x1, x2, period=None, equality_fill_val=0.0, return_pbc_correction=False):
x1 = np.atleast_1d(x1)
x2 = np.atleast_1d(x2)
result = np.sign(x1 - x2)
@@ -956,9 +983,9 @@ def _sign_pbc(x1, x2, period=None, equality_fill_val=0., return_pbc_correction=F
msg = "If period is not None, all values of x and y must be between [0, period)"
raise ValueError(msg)
- d = np.abs(x1-x2)
- pbc_correction = np.sign(period/2. - d)
- result = pbc_correction*result
+ d = np.abs(x1 - x2)
+ pbc_correction = np.sign(period / 2.0 - d)
+ result = pbc_correction * result
if equality_fill_val != 0:
result = np.where(result == 0, equality_fill_val, result)
@@ -970,17 +997,16 @@ def _sign_pbc(x1, x2, period=None, equality_fill_val=0., return_pbc_correction=F
def _relative_positions_and_velocities(x1, x2, period=None, **kwargs):
- s = _sign_pbc(x1, x2, period=period, equality_fill_val=1.)
+ s = _sign_pbc(x1, x2, period=period, equality_fill_val=1.0)
absd = np.abs(x1 - x2)
if period is None:
- xrel = s*absd
+ xrel = s * absd
else:
- xrel = s*np.where(absd > period/2., period - absd, absd)
+ xrel = s * np.where(absd > period / 2.0, period - absd, absd)
try:
- v1 = kwargs['v1']
- v2 = kwargs['v2']
- return xrel, s*(v1-v2)
+ v1 = kwargs["v1"]
+ v2 = kwargs["v2"]
+ return xrel, s * (v1 - v2)
except KeyError:
return xrel
-
diff --git a/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/nfw_profile.py b/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/nfw_profile.py
index d600e7ecd..c547c9c8e 100644
--- a/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/nfw_profile.py
+++ b/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/nfw_profile.py
@@ -20,8 +20,8 @@
from ......sim_manager import sim_defaults
-__all__ = ('NFWProfile', )
-__author__ = ('Andrew Hearin', 'Benedikt Diemer')
+__all__ = ("NFWProfile",)
+__author__ = ("Andrew Hearin", "Benedikt Diemer")
class NFWProfile(AnalyticDensityProf):
@@ -35,13 +35,17 @@ class NFWProfile(AnalyticDensityProf):
implemented in the Halotools code base, see :ref:`nfw_profile_tutorial`.
"""
- def __init__(self,
- cosmology=sim_defaults.default_cosmology,
- redshift=sim_defaults.default_redshift,
- mdef=model_defaults.halo_mass_definition,
- conc_mass_model=model_defaults.conc_mass_model,
- concentration_key=model_defaults.concentration_key,
- **kwargs):
+
+ def __init__(
+ self,
+ cosmology=sim_defaults.default_cosmology,
+ redshift=sim_defaults.default_redshift,
+ mdef=model_defaults.halo_mass_definition,
+ conc_mass_model=model_defaults.conc_mass_model,
+ concentration_key=model_defaults.concentration_key,
+ halo_boundary_key=None,
+ **kwargs
+ ):
r"""
Parameters
----------
@@ -57,6 +61,9 @@ def __init__(self,
String specifying the halo mass definition, e.g., 'vir' or '200m'.
Default is set in `~halotools.empirical_models.model_defaults`.
+ halo_boundary_key : str, optional
+ Default behavior is to use the column associated with the input mdef.
+
conc_mass_model : string or callable, optional
Specifies the function used to model the relation between
NFW concentration and halo mass.
@@ -77,20 +84,24 @@ def __init__(self,
--------
>>> nfw = NFWProfile()
"""
- AnalyticDensityProf.__init__(self, cosmology, redshift, mdef)
+ AnalyticDensityProf.__init__(
+ self, cosmology, redshift, mdef, halo_boundary_key=halo_boundary_key
+ )
- self.gal_prof_param_keys = ['conc_NFWmodel']
- self.halo_prof_param_keys = ['conc_NFWmodel']
+ self.gal_prof_param_keys = ["conc_NFWmodel"]
+ self.halo_prof_param_keys = ["conc_NFWmodel"]
- self.publications = ['arXiv:9611107', 'arXiv:0002395']
+ self.publications = ["arXiv:9611107", "arXiv:0002395"]
- self._initialize_conc_mass_behavior(conc_mass_model,
- concentration_key=concentration_key)
+ self._initialize_conc_mass_behavior(
+ conc_mass_model, concentration_key=concentration_key
+ )
def _initialize_conc_mass_behavior(self, conc_mass_model, **kwargs):
- if conc_mass_model == 'direct_from_halo_catalog':
- self.concentration_key = kwargs.get('concentration_key',
- model_defaults.concentration_key)
+ if conc_mass_model == "direct_from_halo_catalog":
+ self.concentration_key = kwargs.get(
+ "concentration_key", model_defaults.concentration_key
+ )
self.conc_mass_model = conc_mass_model
@@ -170,24 +181,29 @@ def conc_NFWmodel(self, *args, **kwargs):
>>> fake_halo_table = Table({'halo_mvir': fake_masses})
>>> model_conc = nfw.conc_NFWmodel(table=fake_halo_table)
"""
- if self.conc_mass_model == 'direct_from_halo_catalog':
+ if self.conc_mass_model == "direct_from_halo_catalog":
try:
- table = kwargs['table']
+ table = kwargs["table"]
except KeyError:
- msg = ("Must pass ``table`` argument to the ``conc_NFWmodel`` function\n"
- "when ``conc_mass_model`` is set to ``direct_from_halo_catalog``\n")
+ msg = (
+ "Must pass ``table`` argument to the ``conc_NFWmodel`` function\n"
+ "when ``conc_mass_model`` is set to ``direct_from_halo_catalog``\n"
+ )
raise KeyError(msg)
- result = direct_from_halo_catalog(table=table,
- concentration_key=self.concentration_key)
-
- elif self.conc_mass_model == 'dutton_maccio14':
- msg = ("Must either pass a ``prim_haloprop`` argument, \n"
- "or a ``table`` argument with an astropy Table that has the ``{0}`` key")
+ result = direct_from_halo_catalog(
+ table=table, concentration_key=self.concentration_key
+ )
+
+ elif self.conc_mass_model == "dutton_maccio14":
+ msg = (
+ "Must either pass a ``prim_haloprop`` argument, \n"
+ "or a ``table`` argument with an astropy Table that has the ``{0}`` key"
+ )
try:
- mass = kwargs['table'][self.prim_haloprop_key]
+ mass = kwargs["table"][self.prim_haloprop_key]
except:
try:
- mass = kwargs['prim_haloprop']
+ mass = kwargs["prim_haloprop"]
except:
raise KeyError(msg.format(self.prim_haloprop_key))
result = dutton_maccio14(mass, self.redshift)
@@ -473,7 +489,7 @@ def rmax(self, total_mass, conc):
"""
halo_radius = self.halo_mass_to_halo_radius(total_mass)
scale_radius = halo_radius / conc
- return 2.16258*scale_radius
+ return 2.16258 * scale_radius
def vmax(self, total_mass, conc):
r""" Maximum circular velocity of the halo profile,
@@ -608,7 +624,7 @@ def mc_generate_nfw_radial_positions(self, **kwargs):
>>> radial_positions = nfw.mc_generate_nfw_radial_positions(halo_mass = 1e12, conc = 10)
>>> radial_positions = nfw.mc_generate_nfw_radial_positions(halo_radius = 0.25)
"""
- kwargs['mdef'] = self.mdef
- kwargs['cosmology'] = self.cosmology
- kwargs['redshift'] = self.redshift
+ kwargs["mdef"] = self.mdef
+ kwargs["cosmology"] = self.cosmology
+ kwargs["redshift"] = self.redshift
return standalone_mc_generate_nfw_radial_positions(**kwargs)