lenstronomy.LensModel.QuadOptimizer package¶
Submodules¶
lenstronomy.LensModel.QuadOptimizer.multi_plane_fast module¶
- class MultiplaneFast(x_image, y_image, z_lens, z_source, lens_model_list, redshift_list, astropy_instance, param_class, foreground_rays, tol_source=1e-05, numerical_alpha_class=None)[source]¶
Bases:
object
This class accelerates ray tracing computations in multi plane lensing for quadruple image lenses by only computing the deflection from objects in front of the main deflector at z_lens one time.
The first ray tracing computation through the foreground is saved and re-used, but it will always have the same shape as the initial x_image, y_image arrays.
- __init__(x_image, y_image, z_lens, z_source, lens_model_list, redshift_list, astropy_instance, param_class, foreground_rays, tol_source=1e-05, numerical_alpha_class=None)[source]¶
- Parameters:
x_image – x_image to fit
y_image – y_image to fit
z_lens – lens redshift
z_source – source redshift
lens_model_list – list of lens models
redshift_list – list of lens redshifts
astropy_instance – instance of astropy to pass to lens model
param_class – an instance of ParamClass (see documentation in QuadOptimmizer.param_manager)
foreground_rays – (optional) pre-computed foreground rays from a previous iteration, if they are not specified they will be re-computed
tol_source – source plane chi^2 sigma
numerical_alpha_class – class for computing numerically tabulated deflection angles
- chi_square(args_lens, *args, **kwargs)[source]¶
- Parameters:
args_lens – array of lens model parameters being optimized, computed from kwargs_lens in a specified param_class, see documentation in QuadOptimizer.param_manager
- Returns:
total chi^2 penalty (source chi^2 + param chi^2), where param chi^2 is computed by the specified param_class
- logL(args_lens, *args, **kwargs)[source]¶
- Parameters:
args_lens – array of lens model parameters being optimized, computed from kwargs_lens in a specified param_class, see documentation in QuadOptimizer.param_manager
- Returns:
the log likelihood corresponding to the given chi^2
- source_plane_chi_square(args_lens, *args, **kwargs)[source]¶
- Parameters:
args_lens – array of lens model parameters being optimized, computed from kwargs_lens in a specified param_class, see documentation in QuadOptimizer.param_manager
- Returns:
chi2 penalty for the source position (all images must map to the same source coordinate)
- ray_shooting_fast(args_lens)[source]¶
Performs a ray tracing computation through observed coordinates on the sky (self._x_image, self._y_image) to the source plane, returning the final coordinates of each ray on the source plane.
- Parameters:
args_lens – An array of parameters being optimized. The array is computed from a set of key word arguments by an instance of ParamClass (see documentation in QuadOptimizer.param_manager)
- Returns:
the xy coordinate of each ray traced back to the source plane
lenstronomy.LensModel.QuadOptimizer.optimizer module¶
- class Optimizer(x_image, y_image, lens_model_list, redshift_list, z_lens, z_source, parameter_class, astropy_instance=None, numerical_alpha_class=None, particle_swarm=True, re_optimize=False, re_optimize_scale=1.0, pso_convergence_mean=50000, foreground_rays=None, tol_source=1e-05, tol_simplex_func=0.001, simplex_n_iterations=400)[source]¶
Bases:
object
Class which executes the optimization routines. Currently implemented as a particle swarm optimization followed by a downhill simplex routine.
Particle swarm optimizer is modified from the CosmoHammer particle swarm routine with different convergence criteria implemented.
- __init__(x_image, y_image, lens_model_list, redshift_list, z_lens, z_source, parameter_class, astropy_instance=None, numerical_alpha_class=None, particle_swarm=True, re_optimize=False, re_optimize_scale=1.0, pso_convergence_mean=50000, foreground_rays=None, tol_source=1e-05, tol_simplex_func=0.001, simplex_n_iterations=400)[source]¶
- Parameters:
x_image – x_image to fit (should be length 4)
y_image – y_image to fit (should be length 4)
lens_model_list – list of lens models for the system
redshift_list – list of lens redshifts for the system
z_lens – the main deflector redshift, the lens models being optimizer must be at this redshift
z_source – the source redshift
parameter_class – an instance of ParamClass (see documentation in QuadOptimizer.param_manager)
astropy_instance – an instance of astropy to pass to the lens model
numerical_alpha_class – a class to compute numerical deflection angles to pass to the lens model
particle_swarm – bool, whether or not to use a PSO fit first
re_optimize – bool, if True the initial spread of particles will be very tight
re_optimize_scale – float, controls how tight the initial spread of particles is
pso_convergence_mean – when to terminate the PSO fit
foreground_rays – (optional) can pass in pre-computed foreground light rays from a previous fit so as to not waste time recomputing them
tol_source – sigma in the source plane chi^2
tol_simplex_func – tolerance for the downhill simplex optimization
simplex_n_iterations – number of iterations per dimension for the downhill simplex optimization
- optimize(n_particles=50, n_iterations=250, verbose=False, threadCount=1)[source]¶
- Parameters:
n_particles – number of PSO particles, will be ignored if self._particle_swarm is False
n_iterations – number of PSO iterations, will be ignored if self._particle_swarm is False
verbose – whether to print stuff
threadCount – integer; number of threads in multi-threading mode
- Returns:
keyword arguments that map (x_image, y_image) to the same source coordinate (source_x, source_y)
lenstronomy.LensModel.QuadOptimizer.param_manager module¶
- class PowerLawParamManager(kwargs_lens_init)[source]¶
Bases:
object
Base class for handling the translation between key word arguments and parameter arrays for EPL mass models.
This class is intended for use in modeling galaxy-scale lenses
- __init__(kwargs_lens_init)[source]¶
- Parameters:
kwargs_lens_init – the initial kwargs_lens before optimizing
- property to_vary_index¶
The number of lens models being varied in this routine. This is set to 2 because the first three lens models are EPL and SHEAR, and their parameters are being optimized.
The kwargs_list is split at to to_vary_index with indicies < to_vary_index accessed in this class, and lens models with indicies > to_vary_index kept fixed.
Note that this requires a specific ordering of lens_model_list :return:
- bounds(re_optimize, scale=1.0)[source]¶
Sets the low/high parameter bounds for the particle swarm optimization.
NOTE: The low/high values specified here are intended for galaxy-scale lenses. If you want to use this for a different size system you should create a new ParamClass with different settings
- Parameters:
re_optimize – keep a narrow window around each parameter
scale – scales the size of the uncertainty window
- Returns:
- class PowerLawFreeShear(kwargs_lens_init)[source]¶
Bases:
PowerLawParamManager
This class implements a fit of EPL + external shear with every parameter except the power law slope allowed to vary.
- class PowerLawFixedShear(kwargs_lens_init, shear_strength)[source]¶
Bases:
PowerLawParamManager
This class implements a fit of EPL + external shear with every parameter except the power law slope AND the shear strength allowed to vary.
The user should specify shear_strengh in the args_param_class keyword when creating the Optimizer class
- class PowerLawFreeShearMultipole(kwargs_lens_init)[source]¶
Bases:
PowerLawParamManager
This class implements a fit of EPL + external shear + a multipole term with every parameter except the power law slope and multipole moment free to vary.
The mass centroid and orientation of the multipole term are fixed to that of the EPL profile
- property to_vary_index¶
The number of lens models being varied in this routine. This is set to 3 because the first three lens models are EPL, SHEAR, and MULTIPOLE, and their parameters are being optimized.
The kwargs_list is split at to to_vary_index with indicies < to_vary_index accessed in this class, and lens models with indicies > to_vary_index kept fixed.
Note that this requires a specific ordering of lens_model_list :return:
- class PowerLawFixedShearMultipole(kwargs_lens_init, shear_strength)[source]¶
Bases:
PowerLawFixedShear
This class implements a fit of EPL + external shear + a multipole term with every parameter except the power law slope, shear strength, and multipole moment free to vary.
The mass centroid and orientation of the multipole term are fixed to that of the EPL profile
- property to_vary_index¶
The number of lens models being varied in this routine. This is set to 3 because the first three lens models are EPL, SHEAR, and MULTIPOLE, and their parameters are being optimized.
The kwargs_list is split at to to_vary_index with indicies < to_vary_index accessed in this class, and lens models with indicies > to_vary_index kept fixed.
Note that this requires a specific ordering of lens_model_list :return: