lenstronomy.Analysis package¶
Submodules¶
lenstronomy.Analysis.image_reconstruction module¶
- class MultiBandImageReconstruction(multi_band_list, kwargs_model, kwargs_params, multi_band_type='multi-linear', kwargs_likelihood=None, verbose=True)[source]¶
Bases:
object
This class manages the output/results of a fitting process and can conveniently access image reconstruction properties in multi-band fitting. In particular, the fitting result does not come with linear inversion parameters (which may or may not be joint or different for multiple bands) and this class performs the linear inversion for the surface brightness amplitudes and stores them for each individual band to be accessible by the user.
This class is the backbone of the ModelPlot routine that provides the interface of this class with plotting and illustration routines.
- __init__(multi_band_list, kwargs_model, kwargs_params, multi_band_type='multi-linear', kwargs_likelihood=None, verbose=True)[source]¶
- Parameters:
multi_band_list – list of imaging data configuration [[kwargs_data, kwargs_psf, kwargs_numerics], […]]
kwargs_model – model keyword argument list
kwargs_params – keyword arguments of the model parameters, same as output of FittingSequence() ‘kwargs_result’
multi_band_type – string, option when having multiple imaging data sets modelled simultaneously. Options are: - ‘multi-linear’: linear amplitudes are inferred on single data set - ‘linear-joint’: linear amplitudes ae jointly inferred - ‘single-band’: single band
kwargs_likelihood – likelihood keyword arguments as supported by the Likelihood() class
verbose – if True (default), computes and prints the total log-likelihood. This option can be deactivated for speedup purposes (does not run linear inversion again), and reduces the number of prints.
- band_setup(band_index=0)[source]¶
ImageModel() instance and keyword arguments of the model components to execute all the options of the ImSim core modules.
- Parameters:
band_index – integer (>=0) of imaging band in order of multi_band_list input to this class
- Returns:
ImageModel() instance and keyword arguments of the model
- class ModelBand(multi_band_list, kwargs_model, model, error_map, cov_param, param, kwargs_params, image_likelihood_mask_list=None, band_index=0, verbose=True)[source]¶
Bases:
object
Class to plot a single band given the full modeling results This class has its specific role when the linear inference is performed on the joint band level and/or when only a subset of model components get used for this specific band in the modeling.
- __init__(multi_band_list, kwargs_model, model, error_map, cov_param, param, kwargs_params, image_likelihood_mask_list=None, band_index=0, verbose=True)[source]¶
- Parameters:
multi_band_list – list of imaging data configuration [[kwargs_data, kwargs_psf, kwargs_numerics], […]]
kwargs_model – model keyword argument list for the full multi-band modeling
model – 2d numpy array of modeled image for the specified band
error_map – 2d numpy array of size of the image, additional error in the pixels coming from PSF uncertainties
cov_param – covariance matrix of the linear inversion
param – 1d numpy array of the linear coefficients of this imaging band
kwargs_params – keyword argument of keyword argument lists of the different model components selected for the imaging band, NOT including linear amplitudes (not required as being overwritten by the param list)
image_likelihood_mask_list – list of 2d numpy arrays of likelihood masks (for all bands)
band_index – integer of the band to be considered in this class
verbose – if True (default), prints the reduced chi2 value for the current band.
- property model¶
- Returns:
model, 2d numpy array
- property norm_residuals¶
- Returns:
normalized residuals, 2d numpy array
- property image_model_class¶
ImageModel() class instance of the single band with only the model components applied to this band.
- Returns:
SingleBandMultiModel() instance, which inherits the ImageModel instance
- property kwargs_model¶
- Returns:
keyword argument of keyword argument lists of the different model components selected for the imaging band, including linear amplitudes. These format matches the image_model_class() return
- point_source_residuals(aperture_radius)[source]¶
Computes integrated residuals within circular apertures around point sources. This routine can assess the accuracy of point source flux measurements.
- Parameters:
aperture_radius – radius of the aperture considering the residuals around the point sources
- Returns:
list of integrated flux residuals (data - model) within the apertures around the point sources
lenstronomy.Analysis.kinematics_api module¶
- class KinematicsAPI(z_lens, z_source, kwargs_model, kwargs_aperture, kwargs_seeing, anisotropy_model, cosmo=None, lens_model_kinematics_bool=None, light_model_kinematics_bool=None, multi_observations=False, kwargs_numerics_galkin=None, analytic_kinematics=False, Hernquist_approx=False, MGE_light=False, MGE_mass=False, kwargs_mge_light=None, kwargs_mge_mass=None, sampling_number=1000, num_kin_sampling=1000, num_psf_sampling=100)[source]¶
Bases:
object
This class contains routines to compute time delays, magnification ratios, line of sight velocity dispersions etc for a given lens model.
- __init__(z_lens, z_source, kwargs_model, kwargs_aperture, kwargs_seeing, anisotropy_model, cosmo=None, lens_model_kinematics_bool=None, light_model_kinematics_bool=None, multi_observations=False, kwargs_numerics_galkin=None, analytic_kinematics=False, Hernquist_approx=False, MGE_light=False, MGE_mass=False, kwargs_mge_light=None, kwargs_mge_mass=None, sampling_number=1000, num_kin_sampling=1000, num_psf_sampling=100)[source]¶
Initialize the class with the lens model and cosmology.
- Parameters:
z_lens – redshift of lens
z_source – redshift of source
kwargs_model – model keyword arguments, needs ‘lens_model_list’, ‘lens_light_model_list’
kwargs_aperture – spectroscopic aperture keyword arguments, see lenstronomy.Galkin.aperture for options
kwargs_seeing – seeing condition of spectroscopic observation, corresponds to kwargs_psf in the GalKin module specified in lenstronomy.GalKin.psf
cosmo – astropy.cosmology instance, if None then will be set to the default cosmology
lens_model_kinematics_bool – bool list of length of the lens model. Only takes a subset of all the models as part of the kinematics computation ( can be used to ignore substructure, shear etc that do not describe the main deflector potential
light_model_kinematics_bool – bool list of length of the light model. Only takes a subset of all the models as part of the kinematics computation (can be used to ignore light components that do not describe the main deflector
multi_observations – bool, if True uses multi-observation to predict a set of different observations with the GalkinMultiObservation() class. kwargs_aperture and kwargs_seeing require to be lists of the individual observations.
anisotropy_model – type of stellar anisotropy model. See details in MamonLokasAnisotropy() class of lenstronomy.GalKin.anisotropy
analytic_kinematics –
- boolean, if True, used the analytic JAM modeling for
a power-law profile on top of a Hernquist light profile
ATTENTION: This may not be accurate for your specific problem!
Hernquist_approx – bool, if True, uses a Hernquist light profile matched to the half light radius of the deflector light profile to compute the kinematics
MGE_light – bool, if true performs the MGE for the light distribution
MGE_mass – bool, if true performs the MGE for the mass distribution
kwargs_numerics_galkin – numerical settings for the integrated line-of-sight velocity dispersion
kwargs_mge_mass – keyword arguments that go into the MGE decomposition routine
kwargs_mge_light – keyword arguments that go into the MGE decomposition routine
sampling_number – int, number of spectral rendering to compute the light weighted integrated LOS dispersion within the aperture. This keyword should be chosen high enough to result in converged results within the tolerance.
num_kin_sampling – number of kinematic renderings on a total IFU
num_psf_sampling – number of PSF displacements for each kinematic rendering on the IFU
- velocity_dispersion(kwargs_lens, kwargs_lens_light, kwargs_anisotropy, r_eff=None, theta_E=None, gamma=None, kappa_ext=0)[source]¶
API for both, analytic and numerical JAM to compute the velocity dispersion [km/s] This routine uses the galkin_setting() routine for the Galkin configurations (see there what options and input is relevant.
- Parameters:
kwargs_lens – lens model keyword arguments
kwargs_lens_light – lens light model keyword arguments
kwargs_anisotropy – stellar anisotropy keyword arguments
r_eff – projected half-light radius of the stellar light associated with the deflector galaxy, optional, if set to None will be computed in this function with default settings that may not be accurate.
theta_E – Einstein radius (optional)
gamma – power-law slope (optional)
kappa_ext – external convergence (optional)
- Returns:
velocity dispersion [km/s]
- velocity_dispersion_map(kwargs_lens, kwargs_lens_light, kwargs_anisotropy, r_eff=None, theta_E=None, gamma=None, kappa_ext=0, direct_convolve=False, supersampling_factor=1, voronoi_bins=None)[source]¶
API for both, analytic and numerical JAM to compute the velocity dispersion map with IFU data [km/s]
- Parameters:
kwargs_lens – lens model keyword arguments
kwargs_lens_light – lens light model keyword arguments
kwargs_anisotropy – stellar anisotropy keyword arguments
r_eff – projected half-light radius of the stellar light associated with the deflector galaxy, optional, if set to None will be computed in this function with default settings that may not be accurate.
theta_E – circularized Einstein radius, optional, if not provided will either be computed in this function with default settings or not required
gamma – power-law slope at the Einstein radius, optional
kappa_ext – external convergence
direct_convolve – bool, if True, compute the 2D integral numerically
supersampling_factor – supersampling factor for 2D integration grid
voronoi_bins – mapping of the voronoi bins, -1 values for pixels not binned
- Returns:
velocity dispersion map in specified bins or grid in kwargs_aperture, in [km/s] unit
- velocity_dispersion_analytical(theta_E, gamma, r_eff, r_ani, kappa_ext=0)[source]¶
Computes the LOS velocity dispersion of the lens within a slit of size R_slit x dR_slit and seeing psf_fwhm. The assumptions are a Hernquist light profile and the spherical power-law lens model at the first position and an Osipkov and Merritt (‘OM’) stellar anisotropy distribution.
Further information can be found in the AnalyticKinematics() class.
- Parameters:
theta_E – Einstein radius
gamma – power-low slope of the mass profile (=2 corresponds to isothermal)
r_ani – anisotropy radius in units of angles
r_eff – projected half-light radius
kappa_ext – external convergence not accounted in the lens models
- Returns:
velocity dispersion in units [km/s]
- galkin_settings(kwargs_lens, kwargs_lens_light, r_eff=None, theta_E=None, gamma=None)[source]¶
- Parameters:
kwargs_lens – lens model keyword argument list
kwargs_lens_light – deflector light keyword argument list
r_eff – half-light radius (optional)
theta_E – Einstein radius (optional)
gamma – local power-law slope at the Einstein radius (optional)
- Returns:
Galkin() instance and mass and light profiles configured for the Galkin module
- kinematic_lens_profiles(kwargs_lens, MGE_fit=False, model_kinematics_bool=None, theta_E=None, gamma=None, kwargs_mge=None, analytic_kinematics=False)[source]¶
Translates the lenstronomy lens and mass profiles into a (sub) set of profiles that are compatible with the GalKin module to compute the kinematics thereof. The requirement is that the profiles are centered at (0, 0) and that for all profile types there exists a 3d de-projected analytical representation.
- Parameters:
kwargs_lens – lens model parameters
MGE_fit – bool, if true performs the MGE for the mass distribution
model_kinematics_bool – bool list of length of the lens model. Only takes a subset of all the models as part of the kinematics computation (can be used to ignore substructure, shear etc that do not describe the main deflector potential
theta_E – (optional float) estimate of the Einstein radius. If present, does not numerically compute this quantity in this routine numerically
gamma – local power-law slope at the Einstein radius (optional)
kwargs_mge – keyword arguments that go into the MGE decomposition routine
analytic_kinematics – bool, if True, solves the Jeans equation analytically for the power-law mass profile with Hernquist light profile
- Returns:
mass_profile_list, keyword argument list
- kinematic_light_profile(kwargs_lens_light, r_eff=None, MGE_fit=False, model_kinematics_bool=None, Hernquist_approx=False, kwargs_mge=None, analytic_kinematics=False)[source]¶
Setting up of the light profile to compute the kinematics in the GalKin module. The requirement is that the profiles are centered at (0, 0) and that for all profile types there exists a 3d de-projected analytical representation.
- Parameters:
kwargs_lens_light – deflector light model keyword argument list
r_eff – (optional float, else=None) Pre-calculated projected half-light radius of the deflector profile. If not provided, numerical calculation is done in this routine if required.
MGE_fit – boolean, if True performs a Multi-Gaussian expansion of the radial light profile and returns this solution.
model_kinematics_bool – list of booleans to indicate a subset of light profiles to be part of the physical deflector light.
Hernquist_approx – boolean, if True replaces the actual light profile(s) with a Hernquist model with matched half-light radius.
kwargs_mge – keyword arguments that go into the MGE decomposition routine
analytic_kinematics – bool, if True, solves the Jeans equation analytically for the power-law mass profile with Hernquist light profile and adjust the settings accordingly
- Returns:
deflector type list, keyword arguments list
- kinematics_modeling_settings(anisotropy_model, kwargs_numerics_galkin, analytic_kinematics=False, Hernquist_approx=False, MGE_light=False, MGE_mass=False, kwargs_mge_light=None, kwargs_mge_mass=None, sampling_number=1000, num_kin_sampling=1000, num_psf_sampling=100)[source]¶
Return the settings for the kinematic modeling.
- Parameters:
anisotropy_model – type of stellar anisotropy model. See details in MamonLokasAnisotropy() class of lenstronomy.GalKin.anisotropy
analytic_kinematics – boolean, if True, used the analytic JAM modeling for a power-law profile on top of a Hernquist light profile ATTENTION: This may not be accurate for your specific problem!
Hernquist_approx – bool, if True, uses a Hernquist light profile matched to the half light radius of the deflector light profile to compute the kinematics
MGE_light – bool, if true performs the MGE for the light distribution
MGE_mass – bool, if true performs the MGE for the mass distribution
kwargs_numerics_galkin – numerical settings for the integrated line-of-sight velocity dispersion
kwargs_mge_mass – keyword arguments that go into the MGE decomposition routine
kwargs_mge_light – keyword arguments that go into the MGE decomposition routine
sampling_number – number of spectral rendering on a single slit
num_kin_sampling – number of kinematic renderings on a total IFU
num_psf_sampling – number of PSF displacements for each kinematic rendering on the IFU
- Returns:
updated settings
lenstronomy.Analysis.lens_profile module¶
- class LensProfileAnalysis(lens_model)[source]¶
Bases:
object
Class with analysis routines to compute derived properties of the lens model.
- effective_einstein_radius_grid(kwargs_lens, center_x=None, center_y=None, model_bool_list=None, grid_num=200, grid_spacing=0.05, get_precision=False, verbose=True)[source]¶
Computes the radius with mean convergence=1 on a grid.
- Parameters:
kwargs_lens – list of lens model keyword arguments
center_x – position of the center (if not set, is attempting to find it from the parameters kwargs_lens)
center_y – position of the center (if not set, is attempting to find it from the parameters kwargs_lens)
model_bool_list – list of booleans indicating the addition (=True) of a model component in computing the Einstein radius
grid_num – integer, number of grid points to numerically evaluate the convergence and estimate the Einstein radius
grid_spacing – spacing in angular units of the grid
get_precision – If True, return the precision of estimated Einstein radius
verbose (bool) – if True, indicates warning when Einstein radius can not be computed
- Returns:
estimate of the Einstein radius
- effective_einstein_radius(kwargs_lens, r_min=0.001, r_max=10.0, num_points=30)[source]¶
Numerical estimate of the Einstein radius with integral approximation of radial convergence profile.
- Parameters:
kwargs_lens – list of lens model keyword arguments
r_min – minimum radius of the convergence integrand
r_max – maximum radius of the convergence integrand (should be larger than Einstein radius)
num_points – number of radial points in log spacing
- Returns:
estimate of the Einstein radius
- local_lensing_effect(kwargs_lens, ra_pos=0, dec_pos=0, model_list_bool=None)[source]¶
Computes deflection, shear and convergence at (ra_pos,dec_pos) for those part of the lens model not included in the main deflector.
- Parameters:
kwargs_lens – lens model keyword argument list
ra_pos – RA position where to compute the external effect
dec_pos – DEC position where to compute the external effect
model_list_bool – boolean list indicating which models effect to be added to the estimate
- Returns:
alpha_x, alpha_y, kappa, shear1, shear2
- profile_slope(kwargs_lens, radius, center_x=None, center_y=None, model_list_bool=None, num_points=10)[source]¶
Computes the logarithmic power-law slope of a profile. ATTENTION: this is not an observable!
- Parameters:
kwargs_lens – lens model keyword argument list
radius – radius from the center where to compute the logarithmic slope (angular units
center_x – center of profile from where to compute the slope
center_y – center of profile from where to compute the slope
model_list_bool – bool list, indicate which part of the model to consider
num_points – number of estimates around the Einstein radius
- Returns:
logarithmic power-law slope
- mst_invariant_differential(kwargs_lens, radius, center_x=None, center_y=None, model_list_bool=None, num_points=10)[source]¶
Average of the radial stretch differential in radial direction, divided by the radial stretch factor.
\[\xi = \frac{\partial \lambda_{\rm rad}}{\partial r} \frac{1}{\lambda_{\rm rad}}\]This quantity is invariant under the MST. The specific definition is provided by Birrer 2021. Equivalent (proportional) definitions are provided by e.g. Kochanek 2020, Sonnenfeld 2018.
- Parameters:
kwargs_lens – lens model keyword argument list
radius – radius from the center where to compute the MST invariant differential
center_x – center position
center_y – center position
model_list_bool – indicate which part of the model to consider
num_points – number of estimates around the radius
- Returns:
xi
- radial_lens_profile(r_list, kwargs_lens, center_x=None, center_y=None, model_bool_list=None)[source]¶
- Parameters:
r_list – list of radii to compute the spherically averaged lens light profile
center_x – center of the profile
center_y – center of the profile
kwargs_lens – lens parameter keyword argument list
model_bool_list – bool list or None, indicating which profiles to sum over
- Returns:
flux amplitudes at r_list radii azimuthally averaged
- multi_gaussian_lens(kwargs_lens, center_x=None, center_y=None, model_bool_list=None, n_comp=20)[source]¶
Multi-gaussian lens model in convergence space.
- Parameters:
kwargs_lens –
n_comp –
- Returns:
- mass_fraction_within_radius(kwargs_lens, center_x, center_y, theta_E, numPix=100)[source]¶
Computes the mean convergence of all the different lens model components within a spherical aperture.
- Parameters:
kwargs_lens – lens model keyword argument list
center_x – center of the aperture
center_y – center of the aperture
theta_E – radius of aperture
- Returns:
list of average convergences for all the model components
- convergence_peak(kwargs_lens, model_bool_list=None, grid_num=200, grid_spacing=0.01, center_x_init=0, center_y_init=0)[source]¶
Computes the maximal convergence position on a grid and returns its coordinate.
- Parameters:
kwargs_lens – lens model keyword argument list
model_bool_list – bool list (optional) to include certain models or not
- Returns:
center_x, center_y
lenstronomy.Analysis.light2mass module¶
- light2mass_interpol(lens_light_model_list, kwargs_lens_light, numPix=100, deltaPix=0.05, subgrid_res=5, center_x=0, center_y=0)[source]¶
Takes a lens light model and turns it numerically in a lens model (with all lensmodel quantities computed on a grid). Then provides an interpolated grid for the quantities.
- Parameters:
kwargs_lens_light – lens light keyword argument list
numPix – number of pixels per axis for the return interpolation
deltaPix – interpolation/pixel size
center_x – center of the grid
center_y – center of the grid
subgrid_res – subgrid for the numerical integrals
- Returns:
keyword arguments for ‘INTERPOL’ lens model
lenstronomy.Analysis.light_profile module¶
- class LightProfileAnalysis(light_model)[source]¶
Bases:
object
Class with analysis routines to compute derived properties of the lens model.
- ellipticity(kwargs_light, grid_spacing, grid_num, center_x=None, center_y=None, model_bool_list=None, num_iterative=10, iterative=False)[source]¶
Make sure that the window covers all the light, otherwise the moments may give a too low answers.
- Parameters:
kwargs_light – keyword argument list of profiles
center_x – center of profile, if None takes it from the first profile in kwargs_light
center_y – center of profile, if None takes it from the first profile in kwargs_light
model_bool_list – list of booleans to select subsets of the profile
grid_spacing – grid spacing over which the moments are computed
grid_num – grid size over which the moments are computed
iterative (boolean) – if True iteratively adopts an eccentric mask to overcome edge effects
num_iterative (int) – number of iterative changes in ellipticity
- Returns:
eccentricities e1, e2
- half_light_radius(kwargs_light, grid_spacing, grid_num, center_x=None, center_y=None, model_bool_list=None)[source]¶
Computes numerically the half-light-radius of the deflector light and the total photon flux.
- Parameters:
kwargs_light – keyword argument list of profiles
center_x – center of profile, if None takes it from the first profile in kwargs_light
center_y – center of profile, if None takes it from the first profile in kwargs_light
model_bool_list – list of booleans to select subsets of the profile
grid_spacing – grid spacing over which the moments are computed
grid_num – grid size over which the moments are computed
- Returns:
half-light radius
- radial_light_profile(r_list, kwargs_light, center_x=None, center_y=None, model_bool_list=None)[source]¶
- Parameters:
r_list – list of radii to compute the spherically averaged lens light profile
center_x – center of the profile
center_y – center of the profile
kwargs_light – lens light parameter keyword argument list
model_bool_list – bool list or None, indicating which profiles to sum over
- Returns:
flux amplitudes at r_list radii spherically averaged
- multi_gaussian_decomposition(kwargs_light, model_bool_list=None, n_comp=20, center_x=None, center_y=None, r_h=None, grid_spacing=0.02, grid_num=200)[source]¶
Multi-gaussian decomposition of the lens light profile (in 1-dimension)
- Parameters:
kwargs_light – keyword argument list of profiles
center_x – center of profile, if None takes it from the first profile in kwargs_light
center_y – center of profile, if None takes it from the first profile in kwargs_light
model_bool_list – list of booleans to select subsets of the profile
grid_spacing – grid spacing over which the moments are computed for the half-light radius
grid_num – grid size over which the moments are computed
n_comp – maximum number of Gaussian’s in the MGE
r_h – float, half light radius to be used for MGE (optional, otherwise using a numerical grid)
- Returns:
amplitudes, sigmas, center_x, center_y
- multi_gaussian_decomposition_ellipse(kwargs_light, model_bool_list=None, center_x=None, center_y=None, grid_num=100, grid_spacing=0.05, n_comp=20)[source]¶
MGE with ellipticity estimate. Attention: numerical grid settings for ellipticity estimate and radial MGE may not necessarily be the same!
- Parameters:
kwargs_light – keyword argument list of profiles
center_x – center of profile, if None takes it from the first profile in kwargs_light
center_y – center of profile, if None takes it from the first profile in kwargs_light
model_bool_list – list of booleans to select subsets of the profile
grid_spacing – grid spacing over which the moments are computed
grid_num – grid size over which the moments are computed
n_comp – maximum number of Gaussians in the MGE
- Returns:
keyword arguments of the elliptical multi Gaussian profile in lenstronomy conventions
lenstronomy.Analysis.multi_patch_reconstruction module¶
- class MultiPatchReconstruction(multi_band_list, kwargs_model, kwargs_params, multi_band_type='joint-linear', kwargs_likelihood=None, kwargs_pixel_grid=None, verbose=True)[source]¶
Bases:
MultiBandImageReconstruction
This class illustrates the model of disconnected multi-patch modeling with ‘joint-linear’ option in one single array.
- __init__(multi_band_list, kwargs_model, kwargs_params, multi_band_type='joint-linear', kwargs_likelihood=None, kwargs_pixel_grid=None, verbose=True)[source]¶
- Parameters:
multi_band_list – list of imaging data configuration [[kwargs_data, kwargs_psf, kwargs_numerics], […]]
kwargs_model – model keyword argument list
kwargs_params – keyword arguments of the model parameters, same as output of FittingSequence() ‘kwargs_result’
multi_band_type – string, option when having multiple imaging data sets modelled simultaneously. Options are: - ‘multi-linear’: linear amplitudes are inferred on single data set - ‘linear-joint’: linear amplitudes ae jointly inferred - ‘single-band’: single band
kwargs_likelihood – likelihood keyword arguments as supported by the Likelihood() class
kwargs_pixel_grid – keyword argument of PixelGrid() class. This is optional and overwrites a minimal grid Attention for consistent pixel grid definitions!
verbose – if True (default), computes and prints the total log-likelihood. This can de-activated for speedup purposes (does not run linear inversion again), and reduces the number of prints.
- property pixel_grid_joint¶
- Returns:
PixelGrid() class instance covering the entire window of the sky including all individual patches
- image_joint()[source]¶
Patch together the individual patches of data and models.
- Returns:
image_joint, model_joint, norm_residuals_joint
- lens_model_joint()[source]¶
Patch together the individual patches of the lens model (can be discontinues)
- Returns:
2d numpy arrays of kappa_joint, magnification_joint, alpha_x_joint, alpha_y_joint
- source(num_pix, delta_pix, center=None)[source]¶
Source in the same coordinate system as the image.
- Parameters:
num_pix – number of pixels per axes
delta_pix – pixel size
center – list with two entries [center_x, center_y] (optional)
- Returns:
2d surface brightness grid of the reconstructed source and PixelGrid() instance of source grid
lenstronomy.Analysis.td_cosmography module¶
- class TDCosmography(z_lens, z_source, kwargs_model, cosmo_fiducial=None, lens_model_kinematics_bool=None, light_model_kinematics_bool=None, kwargs_seeing=None, kwargs_aperture=None, anisotropy_model=None, **kwargs_kin_api)[source]¶
Bases:
KinematicsAPI
Class equipped to perform a cosmographic analysis from a lens model with added measurements of time delays and kinematics.
This class does not require any cosmological knowledge and can return angular diameter distance estimates self-consistently integrating the kinematics routines and time delay estimates in the lens modeling. This description follows Birrer et al. 2016, 2019.
- __init__(z_lens, z_source, kwargs_model, cosmo_fiducial=None, lens_model_kinematics_bool=None, light_model_kinematics_bool=None, kwargs_seeing=None, kwargs_aperture=None, anisotropy_model=None, **kwargs_kin_api)[source]¶
- Parameters:
z_lens – redshift of deflector
z_source – redshift of source
kwargs_model – model configurations (according to FittingSequence)
cosmo_fiducial – fiducial cosmology used to compute angular diameter distances where required
lens_model_kinematics_bool – (optional) bool list, corresponding to lens models being included into the kinematics modeling
light_model_kinematics_bool – (optional) bool list, corresponding to lens light models being included into the kinematics modeling
kwargs_seeing – seeing conditions (see observation class in Galkin)
kwargs_aperture – aperture keyword arguments (see aperture class in Galkin)
anisotropy_model – string, anisotropy model type
kwargs_kin_api – additional keyword arguments for KinematicsAPI class instance
- time_delays(kwargs_lens, kwargs_ps, kappa_ext=0, original_ps_position=False)[source]¶
Predicts the time delays of the image positions given the fiducial cosmology relative to a straight line without lensing. Negative values correspond to images arriving earlier, and positive signs correspond to images arriving later.
- Parameters:
kwargs_lens – lens model parameters
kwargs_ps – point source parameters
kappa_ext – external convergence (optional)
original_ps_position – boolean (only applies when first point source model is of type ‘LENSED_POSITION’), uses the image positions in the model parameters and does not re-compute images (which might be differently ordered) in case of the lens equation solver
- Returns:
time delays at image positions for the fixed cosmology in units of days
- fermat_potential(kwargs_lens, kwargs_ps, original_ps_position=False)[source]¶
Fermat potential (negative sign means earlier arrival time)
- Parameters:
kwargs_lens – lens model keyword argument list
kwargs_ps – point source keyword argument list
original_ps_position – boolean (only applies when first point source model is of type ‘LENSED_POSITION’), uses the image positions in the model parameters and does not re-compute images (which might be differently ordered) in case of the lens equation solver
- Returns:
Fermat potential of all the image positions in the first point source list entry
- velocity_dispersion_dimension_less(kwargs_lens, kwargs_lens_light, kwargs_anisotropy, r_eff=None, theta_E=None, gamma=None)[source]¶
Sigma**2 = Dd/Dds * c**2 * J(kwargs_lens, kwargs_light, anisotropy) (Equation 4.11 in Birrer et al. 2016 or Equation 6 in Birrer et al. 2019) J() is a dimensionless and cosmological independent quantity only depending on angular units. This function returns J given the lens and light parameters and the anisotropy choice without an external mass sheet correction.
- Parameters:
kwargs_lens – lens model keyword arguments
kwargs_lens_light – lens light model keyword arguments
kwargs_anisotropy – stellar anisotropy keyword arguments
r_eff – projected half-light radius of the stellar light associated with the deflector galaxy, optional, if set to None will be computed in this function with default settings that may not be accurate.
theta_E – pre-computed Einstein radius (optional)
gamma – pre-computed power-law slope of mass profile
- Returns:
dimensionless velocity dispersion (see e.g. Birrer et al. 2016, 2019)
- velocity_dispersion_map_dimension_less(kwargs_lens, kwargs_lens_light, kwargs_anisotropy, r_eff=None, theta_E=None, gamma=None, direct_convolve=False, supersampling_factor=1, voronoi_bins=None)[source]¶
Sigma**2 = Dd/Dds * c**2 * J(kwargs_lens, kwargs_light, anisotropy) (Equation 4.11 in Birrer et al. 2016 or Equation 6 in Birrer et al. 2019) J() is a dimensionless and cosmological independent quantity only depending on angular units. This function returns J given the lens and light parameters and the anisotropy choice without an external mass sheet correction. This routine computes the IFU map of the kinematic quantities.
- Parameters:
kwargs_lens – lens model keyword arguments
kwargs_lens_light – lens light model keyword arguments
kwargs_anisotropy – stellar anisotropy keyword arguments
r_eff – projected half-light radius of the stellar light associated with the deflector galaxy, optional, if set to None will be computed in this function with default settings that may not be accurate.
direct_convolve – bool, if True, compute the 2D integral numerically
supersampling_factor – supersampling factor for 2D integration grid
voronoi_bins – mapping of the voronoi bins, -1 values for pixels not binned
- Returns:
dimensionless velocity dispersion (see e.g. Birrer et al. 2016, 2019)
- static ddt_from_time_delay(d_fermat_model, dt_measured, kappa_s=0, kappa_ds=0, kappa_d=0)[source]¶
Time-delay distance in units of Mpc from the modeled Fermat potential and measured time delay from an image pair.
- Parameters:
d_fermat_model – relative Fermat potential between two images from the same source in units arcsec^2
dt_measured – measured time delay between the same image pair in units of days
kappa_s – external convergence from observer to source
kappa_ds – external convergence from lens to source
kappa_d – external convergence form observer to lens
- Returns:
D_dt, time-delay distance
- static ds_dds_from_kinematics(sigma_v, J, kappa_s=0, kappa_ds=0)[source]¶
Computes the estimate of the ratio of angular diameter distances Ds/Dds from the kinematic estimate of the lens and the measured dispersion.
- Parameters:
sigma_v – velocity dispersion [km/s]
J – dimensionless kinematic constraint (see Birrer et al. 2016, 2019)
- Returns:
Ds/Dds
- ddt_dd_from_time_delay_and_kinematics(d_fermat_model, dt_measured, sigma_v_measured, J, kappa_s=0, kappa_ds=0, kappa_d=0)[source]¶
- Parameters:
d_fermat_model – relative Fermat potential in units arcsec^2
dt_measured – measured relative time delay [days]
sigma_v_measured – 1-sigma Gaussian uncertainty in the measured velocity dispersion
J – modeled dimensionless kinematic estimate
kappa_s – LOS convergence from observer to source
kappa_ds – LOS convergence from deflector to source
kappa_d – LOS convergence from observer to deflector
- Returns:
D_dt, D_d