From 0edb2939df3e677aea1e33bb23515a37edf3a204 Mon Sep 17 00:00:00 2001 From: Andrew Hearin Date: Sun, 2 Feb 2020 09:06:50 -0600 Subject: [PATCH 1/6] Updated installation instructions --- docs/install.rst | 16 ++++++++++------ docs/installing_halotools_with_virtualenv.rst | 6 +++--- 2 files changed, 13 insertions(+), 9 deletions(-) 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:: From ac7db4baca5990890c88e14cdfa15736e094f710 Mon Sep 17 00:00:00 2001 From: Andrew Hearin Date: Fri, 18 Sep 2020 09:07:54 -0500 Subject: [PATCH 2/6] Changes made by atom --- .../hod_models/tests/test_zheng07.py | 71 ++++--- .../factories/hod_mock_factory.py | 184 +++++++++++------- 2 files changed, 152 insertions(+), 103 deletions(-) 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..80d013319 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,28 @@ 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) - - msg = "host_centric_distance key was never mapped during mock population" - assert np.any(model.mock.galaxy_table['host_centric_distance'] > 0), msg + 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 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 From e851791b60de8528c880adee108202d70322ffae Mon Sep 17 00:00:00 2001 From: Andrew Hearin Date: Fri, 18 Sep 2020 09:27:30 -0500 Subject: [PATCH 3/6] Added new unit-test reproducing error reported by @jmsull and @rainwoodman --- .../composite_models/hod_models/tests/test_zheng07.py | 9 +++++++++ 1 file changed, 9 insertions(+) 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 80d013319..79257ad42 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 @@ -104,3 +104,12 @@ def test_host_centric_distance(): 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") + halocat = FakeSim(seed=43) + halocat.halo_table["halo_m200m"] = np.copy(halocat.halo_table["halo_mvir"]) + model.populate_mock(halocat, seed=43) From de727d3073e8a3057b31f705c95b8ce89c027c23 Mon Sep 17 00:00:00 2001 From: Andrew Hearin Date: Fri, 18 Sep 2020 10:10:10 -0500 Subject: [PATCH 4/6] pep8 changes made by atom --- .../centrals/trivial_phase_space.py | 36 +++-- .../analytic_models/profile_model_template.py | 70 +++++--- .../satellites/nfw/biased_nfw_phase_space.py | 149 +++++++++++------- .../satellites/nfw/nfw_profile.py | 76 +++++---- 4 files changed, 202 insertions(+), 129 deletions(-) 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..ba08a8a31 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,13 @@ 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, + **kwargs + ): """ Parameters ---------- @@ -40,11 +42,17 @@ 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 = {} @@ -56,6 +64,6 @@ def __init__(self, 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..49a99c6d5 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) @@ -63,8 +63,8 @@ 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) + cosmology=self.cosmology, redshift=self.redshift, mdef=self.mdef + ) self.halo_boundary_key = model_defaults.get_halo_boundary_key(self.mdef) self.prim_haloprop_key = model_defaults.get_halo_mass_key(self.mdef) @@ -157,11 +157,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 +185,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 +223,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 +271,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 +301,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 +324,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 +359,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 +396,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 +424,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 +452,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 +479,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..a7bcf1561 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 @@ -142,20 +145,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 +176,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 +261,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 +361,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 +389,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 +426,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 +461,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_profile.py b/halotools/empirical_models/phase_space_models/analytic_models/satellites/nfw/nfw_profile.py index d600e7ecd..e7ec7eea0 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,16 @@ 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, + **kwargs + ): r""" Parameters ---------- @@ -79,18 +82,20 @@ def __init__(self, """ AnalyticDensityProf.__init__(self, cosmology, redshift, mdef) - 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 +175,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 +483,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 +618,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) From d3abd0a3ea14e9b8a4152d274c56c80584d071cb Mon Sep 17 00:00:00 2001 From: Andrew Hearin Date: Fri, 18 Sep 2020 10:25:49 -0500 Subject: [PATCH 5/6] Added new keyword argument halo_boundary_key allowing users to override the enforcement of consistent column names between halo mass and halo boundary --- .../composite_models/hod_models/tests/test_zheng07.py | 7 ++++++- .../analytic_models/centrals/trivial_phase_space.py | 6 +++++- .../analytic_models/profile_model_template.py | 10 ++++++++-- .../satellites/nfw/biased_nfw_phase_space.py | 3 +++ .../analytic_models/satellites/nfw/nfw_profile.py | 8 +++++++- 5 files changed, 29 insertions(+), 5 deletions(-) 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 79257ad42..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 @@ -109,7 +109,12 @@ def test_host_centric_distance(): def test_mass_definition_flexibility(): """ Regression test for Issue #993. """ - model = PrebuiltHodModelFactory("zheng07", mdef="200m") + 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/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 ba08a8a31..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 @@ -27,6 +27,7 @@ def __init__( cosmology=sim_defaults.default_cosmology, redshift=sim_defaults.default_redshift, mdef=model_defaults.halo_mass_definition, + halo_boundary_key=None, **kwargs ): """ @@ -59,7 +60,10 @@ def __init__( 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""" 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 49a99c6d5..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 @@ -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 @@ -65,7 +68,10 @@ def __init__(self, cosmology, redshift, mdef, **kwargs): 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) + 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 = [] 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 a7bcf1561..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 @@ -71,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. 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 e7ec7eea0..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 @@ -43,6 +43,7 @@ def __init__( 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""" @@ -60,6 +61,9 @@ def __init__( 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. @@ -80,7 +84,9 @@ def __init__( -------- >>> 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"] From 9856360ac0167b325f6a7ae757df325ea0afbc72 Mon Sep 17 00:00:00 2001 From: Andrew Hearin Date: Fri, 18 Sep 2020 10:26:01 -0500 Subject: [PATCH 6/6] Docstring update --- .../satellites/nfw/nfw_phase_space.py | 112 +++++++++++------- 1 file changed, 69 insertions(+), 43 deletions(-) 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 -