Merge pull request #20 from ajshajib/feature/composite

Hierarchical sampling with composite mass model

hierarc/LensPosterior/base_config.py

@@ -22,6 +22,7 @@ def __init__(
         kwargs_seeing,
         kwargs_numerics_galkin,
         anisotropy_model,
+        lens_model_list=None,
         kwargs_lens_light=None,
         lens_light_model_list=["HERNQUIST"],
         MGE_light=False,
@@ -31,6 +32,7 @@ def __init__(
         num_psf_sampling=100,
         num_kin_sampling=1000,
         multi_observations=False,
+        cosmo_fiducial=None,
     ):
         """
 
@@ -42,27 +44,43 @@ def __init__(
         :param gamma_error: 1-sigma uncertainty on power-law slope
         :param r_eff: half-light radius of the deflector (arc seconds)
         :param r_eff_error: uncertainty on half-light radius
-        :param kwargs_aperture: spectroscopic aperture keyword arguments, see lenstronomy.Galkin.aperture for options
-        :param kwargs_seeing: seeing condition of spectroscopic observation, corresponds to kwargs_psf in the GalKin module specified in lenstronomy.GalKin.psf
-        :param kwargs_numerics_galkin: numerical settings for the integrated line-of-sight velocity dispersion
-        :param anisotropy_model: type of stellar anisotropy model. See details in MamonLokasAnisotropy() class of lenstronomy.GalKin.anisotropy
-        :param multi_observations: bool, if True, interprets kwargs_aperture and kwargs_seeing as lists of multiple
-         observations
+        :param kwargs_aperture: spectroscopic aperture keyword arguments, see
+        lenstronomy.Galkin.aperture for options
+        :param kwargs_seeing: seeing condition of spectroscopic observation, corresponds
+            to kwargs_psf in the GalKin module specified in lenstronomy.GalKin.psf
+        :param kwargs_numerics_galkin: numerical settings for the integrated
+            line-of-sight velocity dispersion
+        :param anisotropy_model: type of stellar anisotropy model. See details in
+            MamonLokasAnisotropy() class of lenstronomy.GalKin.anisotropy
+        :param multi_observations: bool, if True, interprets kwargs_aperture and
+            kwargs_seeing as lists of multiple observations
+        :param lens_model_list: keyword argument list of lens model (optional)
         :param kwargs_lens_light: keyword argument list of lens light model (optional)
-        :param kwargs_mge_light: keyword arguments that go into the MGE decomposition routine
-        :param hernquist_approx: bool, if True, uses the Hernquist approximation for the light profile
+        :param kwargs_mge_light: keyword arguments that go into the MGE decomposition
+            routine
+        :param hernquist_approx: bool, if True, uses the Hernquist approximation for the
+            light profile
+        :param cosmo_fiducial: astropy.cosmology instance, if None,
+            uses astropy's default cosmology
         """
         self._z_lens, self._z_source = z_lens, z_source
-        kwargs_model = {
-            "lens_model_list": ["SPP"],
-            "lens_light_model_list": lens_light_model_list,
-        }
+
+        if lens_model_list is None:
+            kwargs_model = {
+                "lens_model_list": ["SPP"],
+                "lens_light_model_list": lens_light_model_list,
+            }
+        else:
+            kwargs_model = {
+                "lens_model_list": lens_model_list,
+                "lens_light_model_list": lens_light_model_list,
+            }
         TDCosmography.__init__(
             self,
             z_lens,
             z_source,
             kwargs_model,
-            cosmo_fiducial=None,
+            cosmo_fiducial=cosmo_fiducial,
             lens_model_kinematics_bool=None,
             light_model_kinematics_bool=None,
             kwargs_seeing=kwargs_seeing,

hierarc/LensPosterior/kin_constraints.py

@@ -27,13 +27,15 @@ def __init__(
         sigma_v_error_cov_matrix=None,
         kwargs_lens_light=None,
         lens_light_model_list=["HERNQUIST"],
+        lens_model_list=None,
         MGE_light=False,
         kwargs_mge_light=None,
         hernquist_approx=True,
         sampling_number=1000,
         num_psf_sampling=100,
         num_kin_sampling=1000,
         multi_observations=False,
+        cosmo_fiducial=None,
     ):
         """
 
@@ -45,23 +47,36 @@ def __init__(
         :param gamma_error: 1-sigma uncertainty on power-law slope
         :param r_eff: half-light radius of the deflector (arc seconds)
         :param r_eff_error: uncertainty on half-light radius
-        :param sigma_v_measured: numpy array of IFU velocity dispersion of the main deflector in km/s
-        :param sigma_v_error_independent: numpy array of 1-sigma uncertainty in velocity dispersion of the IFU
+        :param sigma_v_measured: numpy array of IFU velocity dispersion of the main
+            deflector in km/s
+        :param sigma_v_error_independent: numpy array of 1-sigma uncertainty in velocity
+            dispersion of the IFU
          observation independent of each other
-        :param sigma_v_error_covariant: covariant error in the measured kinematics shared among all IFU measurements
-        :param sigma_v_error_cov_matrix: error covariance matrix in the sigma_v measurements (km/s)^2
-        :type sigma_v_error_cov_matrix: nxn matrix with n the length of the sigma_v_measured array
-        :param kwargs_aperture: spectroscopic aperture keyword arguments, see lenstronomy.Galkin.aperture for options
-        :param kwargs_seeing: seeing condition of spectroscopic observation, corresponds to kwargs_psf in the GalKin
+        :param sigma_v_error_covariant: covariant error in the measured kinematics
+            shared among all IFU measurements
+        :param sigma_v_error_cov_matrix: error covariance matrix in the sigma_v
+            measurements (km/s)^2
+        :type sigma_v_error_cov_matrix: nxn matrix with n the length of the
+            sigma_v_measured array
+        :param kwargs_aperture: spectroscopic aperture keyword arguments, see
+            lenstronomy.Galkin.aperture for options
+        :param kwargs_seeing: seeing condition of spectroscopic observation, corresponds
+            to kwargs_psf in the GalKin
          module specified in lenstronomy.GalKin.psf
-        :param kwargs_numerics_galkin: numerical settings for the integrated line-of-sight velocity dispersion
-        :param anisotropy_model: type of stellar anisotropy model. See details in MamonLokasAnisotropy() class of
-         lenstronomy.GalKin.anisotropy
+        :param kwargs_numerics_galkin: numerical settings for the integrated
+            line-of-sight velocity dispersion
+        :param anisotropy_model: type of stellar anisotropy model. See details in
+            MamonLokasAnisotropy() class of lenstronomy.GalKin.anisotropy
+        :param lens_model_list: keyword argument list of lens model (optional)
         :param kwargs_lens_light: keyword argument list of lens light model (optional)
-        :param kwargs_mge_light: keyword arguments that go into the MGE decomposition routine
-        :param hernquist_approx: bool, if True, uses the Hernquist approximation for the light profile
-        :param multi_observations: bool, if True, interprets kwargs_aperture and kwargs_seeing as lists of multiple
-         observations
+        :param kwargs_mge_light: keyword arguments that go into the MGE decomposition
+            routine
+        :param hernquist_approx: bool, if True, uses the Hernquist approximation for the
+            light profile
+        :param multi_observations: bool, if True, interprets kwargs_aperture and
+            kwargs_seeing as lists of multiple observations
+        :param cosmo_fiducial: astropy.cosmology instance, if None,
+            uses astropy's default
         """
         self._sigma_v_measured = np.array(sigma_v_measured)
         self._sigma_v_error_independent = np.array(sigma_v_error_independent)
@@ -70,6 +85,7 @@ def __init__(
 
         self._kwargs_lens_light = kwargs_lens_light
         self._anisotropy_model = anisotropy_model
+
         BaseLensConfig.__init__(
             self,
             z_lens,
@@ -84,6 +100,7 @@ def __init__(
             kwargs_seeing,
             kwargs_numerics_galkin,
             anisotropy_model,
+            lens_model_list=lens_model_list,
             kwargs_lens_light=kwargs_lens_light,
             lens_light_model_list=lens_light_model_list,
             MGE_light=MGE_light,
@@ -93,6 +110,7 @@ def __init__(
             num_psf_sampling=num_psf_sampling,
             num_kin_sampling=num_kin_sampling,
             multi_observations=multi_observations,
+            cosmo_fiducial=cosmo_fiducial,
         )
 
     def j_kin_draw(self, kwargs_anisotropy, no_error=False):
@@ -161,9 +179,10 @@ def hierarchy_configuration(self, num_sample_model=20):
     def model_marginalization(self, num_sample_model=20):
         """
 
-        :param num_sample_model: number of samples drawn from the lens and light model posterior to compute the
-         dimensionless kinematic component J()
-        :return: J() as array for each measurement prediction, covariance matrix in sqrt(J)
+        :param num_sample_model: number of samples drawn from the lens and light model
+            posterior to compute the dimensionless kinematic component J()
+        :return: J() as array for each measurement prediction, covariance matrix in
+            sqrt(J)
         """
         num_data = len(self._sigma_v_measured)
         j_kin_matrix = np.zeros(