lenstronomy.ImSim package¶
Subpackages¶
- lenstronomy.ImSim.MultiBand package
- Submodules
- lenstronomy.ImSim.MultiBand.joint_linear module
- lenstronomy.ImSim.MultiBand.multi_data_base module
- lenstronomy.ImSim.MultiBand.multi_linear module
- lenstronomy.ImSim.MultiBand.single_band_multi_model module
SingleBandMultiModel
SingleBandMultiModel.__init__()
SingleBandMultiModel.image()
SingleBandMultiModel.source_surface_brightness()
SingleBandMultiModel.lens_surface_brightness()
SingleBandMultiModel.point_source()
SingleBandMultiModel.image_linear_solve()
SingleBandMultiModel.likelihood_data_given_model()
SingleBandMultiModel.num_param_linear()
SingleBandMultiModel.linear_response_matrix()
SingleBandMultiModel.error_map_source()
SingleBandMultiModel.error_response()
SingleBandMultiModel.update_linear_kwargs()
SingleBandMultiModel.extinction_map()
SingleBandMultiModel.linear_param_from_kwargs()
SingleBandMultiModel.select_kwargs()
- Module contents
- lenstronomy.ImSim.Numerics package
- Submodules
- lenstronomy.ImSim.Numerics.adaptive_numerics module
- lenstronomy.ImSim.Numerics.convolution module
- lenstronomy.ImSim.Numerics.grid module
- lenstronomy.ImSim.Numerics.numba_convolution module
- lenstronomy.ImSim.Numerics.numerics module
- lenstronomy.ImSim.Numerics.partial_image module
- lenstronomy.ImSim.Numerics.point_source_rendering module
- Module contents
Submodules¶
lenstronomy.ImSim.de_lens module¶
- get_param_WLS(A, C_D_inv, d, inv_bool=True)[source]¶
Returns the parameter values given.
- Parameters:
A – response matrix Nd x Ns (Nd = # data points, Ns = # parameters)
C_D_inv – inverse covariance matrix of the data, Nd x Nd, diagonal form
d – data array, 1-d Nd
inv_bool – boolean, whether returning also the inverse matrix or just solve the linear system
- Returns:
1-d array of parameter values
lenstronomy.ImSim.image2source_mapping module¶
- class Image2SourceMapping(lens_model, source_model)[source]¶
Bases:
object
This class handles multiple source planes and performs the computation of predicted surface brightness at given image positions. The class is enable to deal with an arbitrary number of different source planes. There are two different settings:
Single lens plane modelling: In case of a single deflector, lenstronomy models the reduced deflection angles (matched to the source plane in single source plane mode). Each source light model can be added a number (scale_factor) that rescales the reduced deflection angle to the specific source plane.
Multiple lens plane modelling: The multi-plane lens modelling requires the assumption of a cosmology and the redshifts of the multiple lens and source planes. The backwards ray-tracing is performed and stopped at the different source plane redshift to compute the mapping between source to image plane.
- __init__(lens_model, source_model)[source]¶
- Parameters:
lens_model – LensModel() class instance
source_model –
LightModel() class instance. The lightModel includes:
source_scale_factor_list: list of floats corresponding to the rescaled deflection angles to the specific source components. None indicates that the list will be set to 1, meaning a single source plane model (in single lens plane mode).
source_redshift_list: list of redshifts of the light components (in multi lens plane mode)
- property T_ij_start_list¶
List of transverse distances from the observer to the start of the source plane.
- property T_ij_end_list¶
List of transverse distances from the observer to the end of the source plane.
- image2source(x, y, kwargs_lens, index_source, kwargs_special=None)[source]¶
mapping of image plane to source plane coordinates WARNING: for multi lens plane computations and multi source planes, this computation can be slow and should be used as rarely as possible.
- Parameters:
x – image plane coordinate (angle)
y – image plane coordinate (angle)
kwargs_lens – lens model kwargs list
index_source – int, index of source model
- Returns:
source plane coordinate corresponding to the source model of index index_source
- image_flux_joint(x, y, kwargs_lens, kwargs_source, kwargs_special=None, k=None)[source]¶
Computes the surface brightness of all light components at image position (x, y)
- Parameters:
x – coordinate in image plane
y – coordinate in image plane
kwargs_lens – lens model kwargs list
kwargs_source – source model kwargs list
k – None or int or list of int for partial evaluation of light models
- Returns:
surface brightness of all joint light components at image position (x, y)
- image_flux_split(x, y, kwargs_lens, kwargs_source, kwargs_special=None)[source]¶
Computes the surface brightness of all light components at image position (x, y)
- Parameters:
x – coordinate in image plane
y – coordinate in image plane
kwargs_lens – lens model kwargs list
kwargs_source – source model kwargs list
- Returns:
list of responses of every single basis component with default amplitude amp=1, in the same order as the light_model_list
lenstronomy.ImSim.image_linear_solve module¶
- class ImageLinearFit(data_class, psf_class=None, lens_model_class=None, source_model_class=None, lens_light_model_class=None, point_source_class=None, extinction_class=None, kwargs_numerics=None, likelihood_mask=None, psf_error_map_bool_list=None, kwargs_pixelbased=None, linear_solver=True)[source]¶
Bases:
ImageModel
Linear version class, inherits ImageModel.
When light models use pixel-based profile types, such as ‘SLIT_STARLETS’, the WLS linear inversion is replaced by the regularized inversion performed by an external solver. The current pixel-based solver is provided by the SLITronomy plug-in.
- __init__(data_class, psf_class=None, lens_model_class=None, source_model_class=None, lens_light_model_class=None, point_source_class=None, extinction_class=None, kwargs_numerics=None, likelihood_mask=None, psf_error_map_bool_list=None, kwargs_pixelbased=None, linear_solver=True)[source]¶
- Parameters:
data_class – ImageData() instance
psf_class – PSF() instance
lens_model_class – LensModel() instance
source_model_class – LightModel() instance
lens_light_model_class – LightModel() instance
point_source_class – PointSource() instance
kwargs_numerics – keyword arguments passed to the Numerics module
likelihood_mask – 2d boolean array of pixels to be counted in the likelihood calculation/linear optimization
psf_error_map_bool_list – list of boolean of length of point source models. Indicates whether PSF error map is used for the point source model stated as the index.
kwargs_pixelbased – keyword arguments with various settings related to the pixel-based solver (see SLITronomy documentation) being applied to the point sources.
linear_solver – bool, if True (default) fixes the linear amplitude parameters ‘amp’ (avoid sampling) such that they get overwritten by the linear solver solution.
- image_linear_solve(kwargs_lens=None, kwargs_source=None, kwargs_lens_light=None, kwargs_ps=None, kwargs_extinction=None, kwargs_special=None, inv_bool=False)[source]¶
Computes the image (lens and source surface brightness with a given lens model). The linear parameters are computed with a weighted linear least square optimization (i.e. flux normalization of the brightness profiles) However in case of pixel-based modelling, pixel values are constrained by an external solver (e.g. SLITronomy).
- Parameters:
kwargs_lens – list of keyword arguments corresponding to the superposition of different lens profiles
kwargs_source – list of keyword arguments corresponding to the superposition of different source light profiles
kwargs_lens_light – list of keyword arguments corresponding to different lens light surface brightness profiles
kwargs_ps – keyword arguments corresponding to “other” parameters, such as external shear and point source image positions
inv_bool – if True, invert the full linear solver Matrix Ax = y for the purpose of the covariance matrix.
- Returns:
2d array of surface brightness pixels of the optimal solution of the linear parameters to match the data
- image_pixelbased_solve(kwargs_lens=None, kwargs_source=None, kwargs_lens_light=None, kwargs_ps=None, kwargs_extinction=None, kwargs_special=None, init_lens_light_model=None)[source]¶
Computes the image (lens and source surface brightness with a given lens model) using the pixel-based solver.
- Parameters:
kwargs_lens – list of keyword arguments corresponding to the superposition of different lens profiles
kwargs_source – list of keyword arguments corresponding to the superposition of different source light profiles
kwargs_lens_light – list of keyword arguments corresponding to different lens light surface brightness profiles
kwargs_ps – keyword arguments corresponding to point sources
kwargs_extinction – keyword arguments corresponding to dust extinction
kwargs_special – keyword arguments corresponding to “special” parameters
init_lens_light_model – optional initial guess for the lens surface brightness
- Returns:
2d array of surface brightness pixels of the optimal solution of the linear parameters to match the data
- linear_response_matrix(kwargs_lens=None, kwargs_source=None, kwargs_lens_light=None, kwargs_ps=None, kwargs_extinction=None, kwargs_special=None)[source]¶
Computes the linear response matrix (m x n), with n being the data size and m being the coefficients.
- Parameters:
kwargs_lens – lens model keyword argument list
kwargs_source – extended source model keyword argument list
kwargs_lens_light – lens light model keyword argument list
kwargs_ps – point source model keyword argument list
kwargs_extinction – extinction model keyword argument list
kwargs_special – special keyword argument list
- Returns:
linear response matrix
- property data_response¶
Returns the 1d array of the data element that is fitted for (including masking)
- Returns:
1d numpy array
- error_response(kwargs_lens, kwargs_ps, kwargs_special)[source]¶
Returns the 1d array of the error estimate corresponding to the data response.
- Returns:
1d numpy array of response, 2d array of additional errors (e.g. point source uncertainties)
- likelihood_data_given_model(kwargs_lens=None, kwargs_source=None, kwargs_lens_light=None, kwargs_ps=None, kwargs_extinction=None, kwargs_special=None, source_marg=False, linear_prior=None, check_positive_flux=False, linear_solver=True)[source]¶
Computes the likelihood of the data given a model This is specified with the non-linear parameters and a linear inversion and prior marginalisation.
- Parameters:
kwargs_lens – list of keyword arguments corresponding to the superposition of different lens profiles
kwargs_source – list of keyword arguments corresponding to the superposition of different source light profiles
kwargs_lens_light – list of keyword arguments corresponding to different lens light surface brightness profiles
kwargs_ps – keyword arguments corresponding to “other” parameters, such as external shear and point source image positions
kwargs_extinction –
kwargs_special –
source_marg – bool, performs a marginalization over the linear parameters
linear_prior – linear prior width in eigenvalues
check_positive_flux – bool, if True, checks whether the linear inversion resulted in non-negative flux components and applies a punishment in the likelihood if so.
linear_solver – bool, if True (default) fixes the linear amplitude parameters ‘amp’ (avoid sampling) such that they get overwritten by the linear solver solution.
- Returns:
log likelihood (natural logarithm)
- likelihood_data_given_model_solution(model, model_error, cov_matrix, param, kwargs_lens, kwargs_source, kwargs_lens_light, kwargs_ps, source_marg=False, linear_prior=None, check_positive_flux=False)[source]¶
- Parameters:
model –
model_error –
cov_matrix –
param –
kwargs_lens –
kwargs_source –
kwargs_lens_light –
kwargs_ps –
source_marg –
linear_prior –
check_positive_flux –
- Returns:
- num_param_linear(kwargs_lens, kwargs_source, kwargs_lens_light, kwargs_ps)[source]¶
- Returns:
number of linear coefficients to be solved for in the linear inversion
- update_linear_kwargs(param, kwargs_lens, kwargs_source, kwargs_lens_light, kwargs_ps)[source]¶
Links linear parameters to kwargs arguments.
- Parameters:
param – linear parameter vector corresponding to the response matrix
- Returns:
updated list of kwargs with linear parameter values
- linear_param_from_kwargs(kwargs_source, kwargs_lens_light, kwargs_ps)[source]¶
Inverse function of update_linear() returning the linear amplitude list for the keyword argument list.
- Parameters:
kwargs_source –
kwargs_lens_light –
kwargs_ps –
- Returns:
list of linear coefficients
- update_pixel_kwargs(kwargs_source, kwargs_lens_light)[source]¶
Update kwargs arguments for pixel-based profiles with fixed properties such as their number of pixels, scale, and center coordinates (fixed to the origin).
- Parameters:
kwargs_source – list of keyword arguments corresponding to the superposition of different source light profiles
kwargs_lens_light – list of keyword arguments corresponding to the superposition of different lens light profiles
- Returns:
updated kwargs_source and kwargs_lens_light
- reduced_residuals(model, error_map=0)[source]¶
- Parameters:
model – 2d numpy array of the modeled image
error_map – 2d numpy array of additional noise/error terms from model components (such as PSF model uncertainties)
- Returns:
2d numpy array of reduced residuals per pixel
- reduced_chi2(model, error_map=0)[source]¶
Returns reduced chi2 :param model: 2d numpy array of a model predicted image :param error_map: same format as model, additional error component (such as PSF errors) :return: reduced chi2.
- property num_data_evaluate¶
Number of data points to be used in the linear solver :return: number of evaluated data points :rtype: int.
- update_data(data_class)[source]¶
- Parameters:
data_class – instance of Data() class
- Returns:
no return. Class is updated.
- image2array_masked(image)[source]¶
Returns 1d array of values in image that are not masked out for the likelihood computation/linear minimization :param image: 2d numpy array of full image :return: 1d array.
- array_masked2image(array)[source]¶
- Parameters:
array – 1d array of values not masked out (part of linear fitting)
- Returns:
2d array of full image
- error_map_source(kwargs_source, x_grid, y_grid, cov_param)[source]¶
Variance of the linear source reconstruction in the source plane coordinates, computed by the diagonal elements of the covariance matrix of the source reconstruction as a sum of the errors of the basis set.
- Parameters:
kwargs_source – keyword arguments of source model
x_grid – x-axis of positions to compute error map
y_grid – y-axis of positions to compute error map
cov_param – covariance matrix of liner inversion parameters
- Returns:
diagonal covariance errors at the positions (x_grid, y_grid)
- point_source_linear_response_set(kwargs_ps, kwargs_lens, kwargs_special, with_amp=True)[source]¶
- Parameters:
kwargs_ps – point source keyword argument list
kwargs_lens – lens model keyword argument list
kwargs_special – special keyword argument list, may include ‘delta_x_image’ and ‘delta_y_image’
with_amp – bool, if True, relative magnification between multiply imaged point sources are held fixed.
- Returns:
list of positions and amplitudes split in different basis components with applied astrometric corrections
- check_positive_flux(kwargs_source, kwargs_lens_light, kwargs_ps)[source]¶
Checks whether the surface brightness profiles contain positive fluxes and returns bool if True.
- Parameters:
kwargs_source – source surface brightness keyword argument list
kwargs_lens_light – lens surface brightness keyword argument list
kwargs_ps – point source keyword argument list
- Returns:
boolean
lenstronomy.ImSim.image_model module¶
- class ImageModel(data_class, psf_class, lens_model_class=None, source_model_class=None, lens_light_model_class=None, point_source_class=None, extinction_class=None, kwargs_numerics=None, kwargs_pixelbased=None)[source]¶
Bases:
object
This class uses functions of lens_model and source_model to make a lensed image.
- __init__(data_class, psf_class, lens_model_class=None, source_model_class=None, lens_light_model_class=None, point_source_class=None, extinction_class=None, kwargs_numerics=None, kwargs_pixelbased=None)[source]¶
- Parameters:
data_class – instance of ImageData() or PixelGrid() class
psf_class – instance of PSF() class
lens_model_class – instance of LensModel() class
source_model_class – instance of LightModel() class describing the source parameters
lens_light_model_class – instance of LightModel() class describing the lens light parameters
point_source_class – instance of PointSource() class describing the point sources
kwargs_numerics – keyword arguments with various numeric description (see ImageNumerics class for options)
kwargs_pixelbased – keyword arguments with various settings related to the pixel-based solver (see SLITronomy documentation)
- reset_point_source_cache(cache=True)[source]¶
Deletes all the cache in the point source class and saves it from then on.
- Parameters:
cache – boolean, if True, saves the next occuring point source positions in the cache
- Returns:
None
- update_psf(psf_class)[source]¶
Update the instance of the class with a new instance of PSF() with a potentially different point spread function.
- Parameters:
psf_class –
- Returns:
no return. Class is updated.
- source_surface_brightness(kwargs_source, kwargs_lens=None, kwargs_extinction=None, kwargs_special=None, unconvolved=False, de_lensed=False, k=None, update_pixelbased_mapping=True)[source]¶
Computes the source surface brightness distribution.
- Parameters:
kwargs_source – list of keyword arguments corresponding to the superposition of different source light profiles
kwargs_lens – list of keyword arguments corresponding to the superposition of different lens profiles
kwargs_extinction – list of keyword arguments of extinction model
unconvolved – if True: returns the unconvolved light distribution (prefect seeing)
de_lensed – if True: returns the un-lensed source surface brightness profile, otherwise the lensed.
k – integer, if set, will only return the model of the specific index
- Returns:
2d array of surface brightness pixels
- lens_surface_brightness(kwargs_lens_light, unconvolved=False, k=None)[source]¶
Computes the lens surface brightness distribution.
- Parameters:
kwargs_lens_light – list of keyword arguments corresponding to different lens light surface brightness profiles
unconvolved – if True, returns unconvolved surface brightness (perfect seeing), otherwise convolved with PSF kernel
- Returns:
2d array of surface brightness pixels
- point_source(kwargs_ps, kwargs_lens=None, kwargs_special=None, unconvolved=False, k=None)[source]¶
Computes the point source positions and paints PSF convolutions on them.
- Parameters:
kwargs_ps –
kwargs_lens –
kwargs_special –
unconvolved –
k –
- Returns:
- image(kwargs_lens=None, kwargs_source=None, kwargs_lens_light=None, kwargs_ps=None, kwargs_extinction=None, kwargs_special=None, unconvolved=False, source_add=True, lens_light_add=True, point_source_add=True)[source]¶
Make an image with a realisation of linear parameter values “param”.
- Parameters:
kwargs_lens – list of keyword arguments corresponding to the superposition of different lens profiles
kwargs_source – list of keyword arguments corresponding to the superposition of different source light profiles
kwargs_lens_light – list of keyword arguments corresponding to different lens light surface brightness profiles
kwargs_ps – keyword arguments corresponding to “other” parameters, such as external shear and point source image positions
unconvolved – if True: returns the unconvolved light distribution (prefect seeing)
source_add – if True, compute source, otherwise without
lens_light_add – if True, compute lens light, otherwise without
point_source_add – if True, add point sources, otherwise without
- Returns:
2d array of surface brightness pixels of the simulation
- extinction_map(kwargs_extinction=None, kwargs_special=None)[source]¶
Differential extinction per pixel.
- Parameters:
kwargs_extinction – list of keyword arguments corresponding to the optical depth models tau, such that extinction is exp(-tau)
kwargs_special – keyword arguments, additional parameter to the extinction
- Returns:
2d array of size of the image