Merge branch 'main' into feature/array_testing

Author: rhayes777

autolens/__init__.py

@@ -124,4 +124,4 @@
 
 conf.instance.register(__file__)
 
-__version__ = "2024.11.6.1"
+__version__ = "2024.11.13.2"

autolens/aggregator/tracer.py

@@ -58,17 +58,18 @@ def _tracer_from(
 
     tracer = Tracer(galaxies=galaxies, cosmology=cosmology)
 
-    if len(fit.children) > 0:
-        logger.info(
-            """
-            Using database for a fit with multiple summed Analysis objects.
-
-            Tracer objects do not fully support this yet (e.g. model parameters which vary over analyses may be incorrect)
-            so proceed with caution!
-            """
-        )
-
-        return [tracer] * len(fit.children)
+    if fit.children is not None:
+        if len(fit.children) > 0:
+            logger.info(
+                """
+                Using database for a fit with multiple summed Analysis objects.
+    
+                Tracer objects do not fully support this yet (e.g. model parameters which vary over analyses may be incorrect)
+                so proceed with caution!
+                """
+            )
+
+            return [tracer] * len(fit.children)
 
     return [tracer]
 

autolens/analysis/result.py

@@ -1,3 +1,4 @@
+import logging
 import os
 import numpy as np
 from typing import Optional, Union
@@ -15,6 +16,8 @@
 from autolens.lens.tracer import Tracer
 from autolens.point.solver import PointSolver
 
+logger = logging.getLogger(__name__)
+
 
 class Result(AgResultDataset):
     @property
@@ -101,8 +104,76 @@ def image_plane_multiple_image_positions(self) -> aa.Grid2DIrregular:
             source_plane_coordinate=self.source_plane_centre.in_list[0],
         )
 
+        if multiple_images.shape[0] == 1:
+            return self.image_plane_multiple_image_positions_for_single_image_from()
+
         return aa.Grid2DIrregular(values=multiple_images)
 
+    def image_plane_multiple_image_positions_for_single_image_from(
+        self, increments: int = 20
+    ) -> aa.Grid2DIrregular:
+        """
+        If the standard point solver only locates one multiple image, finds one or more additional images, which are
+        not technically multiple image in the point source regime, but are close enough to it they can be used
+        in a position threshold likelihood.
+
+        This is performed by incrementally moving the source-plane centre's coordinates towards the centre of the
+        source-plane at (0.0", 0.0"). This ensures that the centre will eventually go inside the caustic, where
+        multiple images are formed.
+
+        To move the source-plane centre, the original source-plane centre is multiplied by a factor that decreases
+        from 1.0 to 0.0 in increments of 1/increments. For example, if the source-plane centre is (1.0", -0.5") and
+        the `factor` is 0.5, the input source-plane centre is (0.5", -0.25").
+
+        The multiple images are always computed for the same mass model, thus they will always be valid multiple images
+        for the model being fitted, but as the factor decrease the multiple images may move furhter from their observed
+        positions.
+
+        Parameters
+        ----------
+        increments
+            The number of increments the source-plane centre is moved to compute multiple images.
+        """
+
+        logger.info(
+            """
+        Could not find multiple images for maximum likelihood lens model.
+        
+        Incrementally moving source centre inwards towards centre of source-plane until caustic crossing occurs
+        and multiple images are formed.        
+        """
+        )
+
+        grid = self.analysis.dataset.mask.derive_grid.all_false
+
+        centre = self.source_plane_centre.in_list[0]
+
+        solver = PointSolver.for_grid(
+            grid=grid,
+            pixel_scale_precision=0.001,
+        )
+
+        for i in range(1, increments):
+            factor = 1.0 - (1.0 * (i / increments))
+
+            multiple_images = solver.solve(
+                tracer=self.max_log_likelihood_tracer,
+                source_plane_coordinate=(centre[0] * factor, centre[1] * factor),
+            )
+
+            if multiple_images.shape[0] > 1:
+                return aa.Grid2DIrregular(values=multiple_images)
+
+        logger.info(
+            """
+        Could not find multiple images for maximum likelihood lens model, even after incrementally moving the source
+        centre inwards to the centre of the source-plane.
+
+        Set the multiple image postiions to two images at (1.0", 1.0") so code continues to run.
+        """
+        )
+        return aa.Grid2DIrregular(values=[(1.0, 1.0), (1.0, 1.0)])
+
     def positions_threshold_from(
         self,
         factor=1.0,
@@ -168,6 +239,7 @@ def positions_likelihood_from(
         minimum_threshold=None,
         use_resample=False,
         positions: Optional[aa.Grid2DIrregular] = None,
+        mass_centre_radial_distance_min: float = None,
     ) -> Union[PositionsLHPenalty, PositionsLHResample]:
         """
         Returns a `PositionsLH` object from the result of a lens model-fit, where the maximum log likelihood mass
@@ -178,6 +250,11 @@ def positions_likelihood_from(
         parametric source) can be used to determine the multiple image positions and threshold for a more complex
         subsequent fit (e.g. power-law mass model, pixelized source).
 
+        The mass model central image is removed from the solution, as this is rarely physically observed and therefore
+        should not be included in the likelihood penalty or resampling. It is removed by setting a positive
+        magnification threshold in the `PointSolver`. For strange lens models the central image may still be
+        solved for, in which case the `mass_centre_radial_distance_min` parameter can be used to remove it.
+
         Parameters
         ----------
         factor
@@ -193,6 +270,10 @@ def positions_likelihood_from(
         positions
             If input, these positions are used instead of the computed multiple image positions from the lens mass
             model.
+        mass_centre_radial_distance_min
+            The minimum radial distance from the mass model centre that a multiple image position must be to be
+            included in the likelihood penalty or resampling. If `None` all positions are used. This is an additional
+            method to remove central images that may make it through the point solver's magnification threshold.
 
         Returns
         -------
@@ -207,6 +288,18 @@ def positions_likelihood_from(
             if positions is None
             else positions
         )
+
+        if mass_centre_radial_distance_min is not None:
+            mass_centre = self.max_log_likelihood_tracer.extract_attribute(
+                cls=ag.mp.MassProfile, attr_name="centre"
+            )
+
+            distances = positions.distances_to_coordinate_from(
+                coordinate=mass_centre[0]
+            )
+
+            positions = positions[distances > mass_centre_radial_distance_min]
+
         threshold = self.positions_threshold_from(
             factor=factor, minimum_threshold=minimum_threshold, positions=positions
         )

autolens/config/visualize/plots.yaml

@@ -55,7 +55,7 @@
     all_at_end_png: true                     # Plot all individual plots listed below as .png (even if False)?
     all_at_end_fits: true                    # Plot all individual plots listed below as .fits (even if False)?
     all_at_end_pdf: false                    # Plot all individual plots listed below as publication-quality .pdf (even if False)?
-    errors: false
+    reconstruction_noise_map: false
     reconstructed_image: false
     reconstruction: false
     regularization_weights: false

autolens/imaging/plot/fit_imaging_plotters.py

@@ -148,7 +148,7 @@ def figures_2d_of_planes(
         subtracted_image: bool = False,
         model_image: bool = False,
         plane_image: bool = False,
-        plane_errors: bool = False,
+        plane_noise_map: bool = False,
         plane_signal_to_noise_map: bool = False,
         use_source_vmax: bool = False,
         zoom_to_brightest: bool = True,
@@ -180,12 +180,12 @@ def figures_2d_of_planes(
             Whether to make a 2D plot (via `imshow`) of the image of a plane in its source-plane (e.g. unlensed).
             Depending on how the fit is performed, this could either be an image of light profiles of the reconstruction
             of an `Inversion`.
-        plane_errors
-            Whether to make a 2D plot (via `imshow`) of the errors of a plane in its source-plane, where the
-            errors can only be computed when a pixelized source reconstruction is performed and they correspond to
-            the errors in each reconstructed pixel as given by the inverse curvature matrix.
+        plane_noise_map
+            Whether to make a 2D plot of the noise-map of a plane in its source-plane, where the
+            noise map can only be computed when a pixelized source reconstruction is performed and they correspond to
+            the noise map in each reconstructed pixel as given by the inverse curvature matrix.
         plane_signal_to_noise_map
-            Whether to make a 2D plot (via `imshow`) of the signal-to-noise map of a plane in its source-plane,
+            Whether to make a 2D plot of the signal-to-noise map of a plane in its source-plane,
             where the signal-to-noise map values can only be computed when a pixelized source reconstruction and they
             are the ratio of reconstructed flux to error in each pixel.
         use_source_vmax
@@ -295,12 +295,6 @@ def figures_2d_of_planes(
 
                     pix = self.tracer.planes[plane_index].cls_list_from(cls=aa.Pixelization)[0]
 
-                    if isinstance(pix.mesh, aa.mesh.Delaunay):
-                        try:
-                            self.mat_plot_2d.cmap.kwargs.pop("vmax")
-                        except KeyError:
-                            pass
-
                     inversion_plotter = self.inversion_plotter_of_plane(
                         plane_index=plane_index
                     )
@@ -318,7 +312,7 @@ def figures_2d_of_planes(
                 except KeyError:
                     pass
 
-            if plane_errors:
+            if plane_noise_map:
 
                 if self.tracer.planes[plane_index].has(cls=aa.Pixelization):
 
@@ -328,7 +322,7 @@ def figures_2d_of_planes(
 
                     inversion_plotter.figures_2d_of_pixelization(
                         pixelization_index=0,
-                        errors=True,
+                        reconstruction_noise_map=True,
                         zoom_to_brightest=zoom_to_brightest,
                         interpolate_to_uniform=interpolate_to_uniform
                     )

autolens/imaging/simulator.py

@@ -9,21 +9,30 @@ class SimulatorImaging(aa.SimulatorImaging):
 
     def via_tracer_from(self, tracer : Tracer, grid : aa.type.Grid2DLike) -> aa.Imaging:
         """
-        Simulate an `Imaging` dataset from an input tracer and grid.
+        Simulate an `Imaging` dataset from an input `Tracer` object and a 2D grid of (y,x) coordinates.
 
-        The tracer is used to perform ray-tracing and generate the image of the strong lens galaxies (e.g.
-        the lens light, lensed source light, etc) which is simulated.
+        The mass profiles of each galaxy in the tracer are used to perform ray-tracing of the input 2D grid and
+        their light profiles are used to generate the image of the galaxies which are simulated.
 
         The steps of the `SimulatorImaging` simulation process (e.g. PSF convolution, noise addition) are
-        described in the `SimulatorImaging` `__init__` method docstring.
+        described in the `SimulatorImaging` `__init__` method docstring, found in the PyAutoArray project.
+
+        If one of more galaxy light profiles are a `LightProfileSNR` object, the `intensity` of the light profile is
+        automatically set such that the signal-to-noise ratio of the light profile is equal to its input
+        `signal_to_noise_ratio` value.
+
+        For example, if a `LightProfileSNR` object has a `signal_to_noise_ratio` of 5.0, the intensity of the light
+        profile is set such that the peak surface brightness of the profile is 5.0 times the background noise level of
+        the image.
 
         Parameters
         ----------
         tracer
             The tracer, which describes the ray-tracing and strong lens configuration used to simulate the imaging
-            dataset.
+            dataset as well as the light profiles of the galaxies used to simulate the image of the galaxies.
         grid
-            The image-plane grid which the image of the strong lens is generated on.
+            The 2D grid of (y,x) coordinates which the mass profiles of the galaxies in the tracer are ray-traced using
+            in order to generate the image of the galaxies via their light profiles.
         """
 
         tracer.set_snr_of_snr_light_profiles(
@@ -44,21 +53,30 @@ def via_tracer_from(self, tracer : Tracer, grid : aa.type.Grid2DLike) -> aa.Imag
 
     def via_galaxies_from(self, galaxies : List[ag.Galaxy], grid : aa.type.Grid2DLike) -> aa.Imaging:
         """
-        Simulate an `Imaging` dataset from an input list of galaxies and grid.
+        Simulate an `Imaging` dataset from an input list of `Galaxy` objects and a 2D grid of (y,x) coordinates.
 
-        The galaxies are used to create a tracer, which performs ray-tracing and generate the image of the strong lens
-        galaxies (e.g. the lens light, lensed source light, etc) which is simulated.
+        The galaxies are used to create a `Tracer`. The mass profiles of each galaxy in the tracer are used to
+        perform ray-tracing of the input 2D grid and their light profiles are used to generate the image of the
+        galaxies which are simulated.
 
         The steps of the `SimulatorImaging` simulation process (e.g. PSF convolution, noise addition) are
         described in the `SimulatorImaging` `__init__` method docstring.
 
+        If one of more galaxy light profiles are a `LightProfileSNR` object, the `intensity` of the light profile is
+        automatically set such that the signal-to-noise ratio of the light profile is equal to its input
+        `signal_to_noise_ratio` value.
+
+        For example, if a `LightProfileSNR` object has a `signal_to_noise_ratio` of 5.0, the intensity of the light
+        profile is set such that the peak surface brightness of the profile is 5.0 times the background noise level of
+        the image.
+
         Parameters
         ----------
         galaxies
             The galaxies used to create the tracer, which describes the ray-tracing and strong lens configuration
             used to simulate the imaging dataset.
         grid
-            The image-plane grid which the image of the strong lens is generated on.
+            The image-plane 2D grid of (y,x) coordinates grid which the image of the strong lens is generated on.
         """
 
         tracer = Tracer(galaxies=galaxies)
@@ -87,7 +105,7 @@ def via_deflections_and_galaxies_from(self, deflections : aa.VectorYX2D, galaxie
             The galaxies used to create the tracer, which describes the ray-tracing and strong lens configuration
             used to simulate the imaging dataset.
         grid
-            The image-plane grid which the image of the strong lens is generated on.
+            The image-plane 2D grid of (y,x) coordinates grid which the image of the strong lens is generated on.
         """
         grid = aa.Grid2D.uniform(
             shape_native=deflections.shape_native,
@@ -105,10 +123,10 @@ def via_source_image_from(self, tracer : Tracer, grid : aa.type.Grid2DLike, sour
         Simulate an `Imaging` dataset from an input image of a source galaxy.
 
         This input image is on a uniform and regular 2D array, meaning it can simulate the source's irregular
-        and assymetric source galaxy morphological features.
+        and asymmetric source galaxy morphological features.
 
         The typical use case is inputting the image of an irregular galaxy in the source-plane (whose values are
-        on a uniform array) and using this function computing the lensed image of this source galaxy.
+        on a uniform array) and using this function to compute the lensed image of this source galaxy.
 
         The tracer is used to perform ray-tracing and generate the image of the strong lens galaxies (e.g.
         the lens light, lensed source light, etc) which is simulated.
@@ -125,7 +143,7 @@ def via_source_image_from(self, tracer : Tracer, grid : aa.type.Grid2DLike, sour
             The tracer, which describes the ray-tracing and strong lens configuration used to simulate the imaging
             dataset.
         grid
-            The image-plane grid which the image of the strong lens is generated on.
+            The image-plane 2D grid of (y,x) coordinates grid which the image of the strong lens is generated on.
         source_image
             The image of the source-plane and source galaxy which is interpolated to compute the lensed image.
         """

autolens/interferometer/plot/fit_interferometer_plotters.py

@@ -375,7 +375,7 @@ def figures_2d_of_planes(
         self,
         plane_index: Optional[int] = None,
         plane_image: bool = False,
-        plane_errors: bool = False,
+        plane_noise_map: bool = False,
         plane_signal_to_noise_map: bool = False,
         zoom_to_brightest: bool = True,
         interpolate_to_uniform: bool = False,
@@ -398,12 +398,12 @@ def figures_2d_of_planes(
             Whether to make a 2D plot (via `imshow`) of the image of a plane in its source-plane (e.g. unlensed).
             Depending on how the fit is performed, this could either be an image of light profiles of the reconstruction
             of an `Inversion`.
-        plane_errors
-            Whether to make a 2D plot (via `imshow`) of the errors of a plane in its source-plane, where the
-            errors can only be computed when a pixelized source reconstruction is performed and they correspond to
-            the errors in each reconstructed pixel as given by the inverse curvature matrix.
+        plane_noise_map
+            Whether to make a 2D plot of the noise-map of a plane in its source-plane, where the
+            noise map can only be computed when a pixelized source reconstruction is performed and they correspond to
+            the noise map in each reconstructed pixel as given by the inverse curvature matrix.
         plane_signal_to_noise_map
-            Whether to make a 2D plot (via `imshow`) of the signal-to-noise map of a plane in its source-plane,
+            Whether to make a 2D plot of the signal-to-noise map of a plane in its source-plane,
             where the signal-to-noise map values can only be computed when a pixelized source reconstruction and they
             are the ratio of reconstructed flux to error in each pixel.
         zoom_to_brightest
@@ -430,15 +430,15 @@ def figures_2d_of_planes(
                     interpolate_to_uniform=interpolate_to_uniform,
                 )
 
-        if plane_errors:
+        if plane_noise_map:
             if self.tracer.planes[plane_index].has(cls=aa.Pixelization):
                 inversion_plotter = self.inversion_plotter_of_plane(
                     plane_index=plane_index
                 )
 
                 inversion_plotter.figures_2d_of_pixelization(
                     pixelization_index=0,
-                    errors=True,
+                    reconstruction_noise_map=True,
                     zoom_to_brightest=zoom_to_brightest,
                     interpolate_to_uniform=interpolate_to_uniform,
                 )