Merge pull request #257 from Jammy2211/feature/sky_refactor

Sky Background Modeling

Author: Jammy2211

autolens/imaging/fit_imaging.py

@@ -21,6 +21,7 @@ def __init__(
         self,
         dataset: aa.Imaging,
         tracer: Tracer,
+        sky: Optional[ag.LightProfile] = None,
         adapt_images: Optional[ag.AdaptImages] = None,
         settings_inversion: aa.SettingsInversion = aa.SettingsInversion(),
         preloads: Preloads = Preloads(),
@@ -56,6 +57,8 @@ def __init__(
             The imaging dataset which is fitted by the galaxies in the tracer.
         tracer
             The tracer of galaxies whose light profile images are used to fit the imaging data.
+        sky
+            Model component used to represent the background sky emission in an image (e.g. a `Sky` light profile).
         adapt_images
             Contains the adapt-images which are used to make a pixelization's mesh and regularization adapt to the
             reconstructed galaxy's morphology.
@@ -71,9 +74,10 @@ def __init__(
 
         super().__init__(dataset=dataset, run_time_dict=run_time_dict)
         AbstractFitInversion.__init__(
-            self=self, model_obj=tracer, settings_inversion=settings_inversion
+            self=self, model_obj=tracer, sky=sky, settings_inversion=settings_inversion
         )
 
+        self.sky = sky
         self.tracer = tracer
 
         self.adapt_images = adapt_images
@@ -92,11 +96,17 @@ def blurred_image(self) -> aa.Array2D:
 
         if self.preloads.blurred_image is None:
 
+            if isinstance(self.sky, ag.lp.Sky):
+                image = self.sky.image_2d_from(grid=self.dataset.grid)
+            else:
+                image = np.zeros(self.dataset.shape_slim)
+
             return self.tracer.blurred_image_2d_from(
                 grid=self.dataset.grid,
                 convolver=self.dataset.convolver,
                 blurring_grid=self.dataset.blurring_grid,
-            )
+            ) + image
+
         return self.preloads.blurred_image
 
     @property
@@ -111,6 +121,7 @@ def tracer_to_inversion(self) -> TracerToInversion:
 
         return TracerToInversion(
             tracer=self.tracer,
+            sky=self.sky,
             dataset=self.dataset,
             data=self.profile_subtracted_image,
             noise_map=self.noise_map,

autolens/imaging/model/analysis.py

@@ -156,13 +156,16 @@ def fit_from(
             instance=instance, run_time_dict=run_time_dict
         )
 
+        sky = self.sky_via_instance_from(instance=instance)
+
         adapt_images = self.adapt_images_via_instance_from(instance=instance)
 
         preloads = preload_overwrite or self.preloads
 
         return FitImaging(
             dataset=self.dataset,
             tracer=tracer,
+            sky=sky,
             adapt_images=adapt_images,
             settings_inversion=self.settings_inversion,
             preloads=preloads,

autolens/imaging/plot/fit_imaging_plotters.py

@@ -467,7 +467,11 @@ def subplot_fit_log10(self):
 
         self.set_title(label=None)
 
-        self.figures_2d(signal_to_noise_map=True)
+        try:
+            self.figures_2d(signal_to_noise_map=True)
+        except ValueError:
+            pass
+
         self.figures_2d(model_image=True)
 
         self.set_title(label="Lens Light Model Image")

autolens/interferometer/fit_interferometer.py

@@ -84,7 +84,7 @@ def __init__(
 
         super().__init__(dataset=dataset, run_time_dict=run_time_dict)
         AbstractFitInversion.__init__(
-            self=self, model_obj=tracer, settings_inversion=settings_inversion
+            self=self, model_obj=tracer, sky=None, settings_inversion=settings_inversion
         )
 
     @property

autolens/interferometer/model/analysis.py

@@ -328,11 +328,7 @@ def visualize(self, paths: af.DirectoryPaths, instance, during_analysis):
             except IndexError:
                 pass
 
-    def make_result(
-        self,
-        samples: af.SamplesPDF,
-            search_internal = None
-    ):
+    def make_result(self, samples: af.SamplesPDF, search_internal=None):
         """
         After the non-linear search is complete create its `Result`, which includes:
 
@@ -358,7 +354,9 @@ def make_result(
         ResultImaging
             The result of fitting the model to the imaging dataset, via a non-linear search.
         """
-        return ResultInterferometer(samples=samples, analysis=self, search_internal=search_internal)
+        return ResultInterferometer(
+            samples=samples, analysis=self, search_internal=search_internal
+        )
 
     def save_attributes(self, paths: af.DirectoryPaths):
         """

autolens/lens/to_inversion.py

@@ -6,7 +6,7 @@
 import autogalaxy as ag
 
 from autoarray.inversion.inversion.factory import inversion_unpacked_from
-
+from autogalaxy.profiles.light.basis import Basis
 from autolens.analysis.preloads import Preloads
 
 
@@ -18,6 +18,7 @@ def __init__(
         data: Optional[Union[aa.Array2D, aa.Visibilities]] = None,
         noise_map: Optional[Union[aa.Array2D, aa.VisibilitiesNoiseMap]] = None,
         w_tilde: Optional[Union[aa.WTildeImaging, aa.WTildeInterferometer]] = None,
+        sky: Optional[Basis] = None,
         adapt_images: Optional[ag.AdaptImages] = None,
         settings_inversion: aa.SettingsInversion = aa.SettingsInversion(),
         preloads=Preloads(),
@@ -30,6 +31,7 @@ def __init__(
             data=data,
             noise_map=noise_map,
             w_tilde=w_tilde,
+            sky=sky,
             adapt_images=adapt_images,
             settings_inversion=settings_inversion,
             preloads=preloads,
@@ -70,6 +72,7 @@ def lp_linear_func_list_galaxy_dict(
         for plane_index, galaxies in enumerate(self.planes):
             plane_to_inversion = ag.GalaxiesToInversion(
                 galaxies=galaxies,
+                sky=self.sky,
                 dataset=self.dataset,
                 grid=traced_grids_of_planes_list[plane_index],
                 blurring_grid=traced_blurring_grids_of_planes_list[plane_index],
@@ -133,6 +136,7 @@ def image_plane_mesh_grid_pg_list(self) -> List[List]:
         for galaxies in self.planes:
             to_inversion = ag.GalaxiesToInversion(
                 galaxies=galaxies,
+                sky=self.sky,
                 grid_pixelization=self.dataset.grid,
                 noise_map=self.noise_map,
                 adapt_images=self.adapt_images,
@@ -206,6 +210,7 @@ def mapper_galaxy_dict(self) -> Dict[aa.AbstractMapper, ag.Galaxy]:
             if galaxies.has(cls=aa.Pixelization):
                 to_inversion = ag.GalaxiesToInversion(
                     galaxies=galaxies,
+                    sky=self.sky,
                     grid_pixelization=traced_grids_of_planes_list[plane_index],
                     preloads=self.preloads,
                     noise_map=self.noise_map,

autolens/point/model/analysis.py

@@ -101,9 +101,11 @@ def visualize(self, paths, instance, during_analysis):
     def make_result(
         self,
         samples: af.SamplesPDF,
-            search_internal = None,
+        search_internal=None,
     ):
-        return ResultPoint(samples=samples, analysis=self, search_internal=search_internal)
+        return ResultPoint(
+            samples=samples, analysis=self, search_internal=search_internal
+        )
 
     def save_attributes(self, paths: af.DirectoryPaths):
         self.point_dict.output_to_json(

autolens/quantity/model/analysis.py

@@ -127,7 +127,7 @@ def visualize(
     def make_result(
         self,
         samples: af.SamplesPDF,
-        search_internal = None,
+        search_internal=None,
     ) -> ResultQuantity:
         """
         After the non-linear search is complete create its `ResultQuantity`, which includes:
@@ -156,4 +156,6 @@ def make_result(
         ResultQuantity
             The result of fitting the model to the imaging dataset, via a non-linear search.
         """
-        return ResultQuantity(samples=samples, analysis=self, search_internal=search_internal)
+        return ResultQuantity(
+            samples=samples, analysis=self, search_internal=search_internal
+        )

test_autolens/imaging/test_fit_imaging.py

@@ -1,3 +1,4 @@
+import copy
 import numpy as np
 import pytest
 
@@ -231,6 +232,60 @@ def test__fit_figure_of_merit(masked_imaging_7x7, masked_imaging_covariance_7x7)
     assert fit.figure_of_merit == pytest.approx(-341415.8258823, 1.0e-4)
 
 
+def test__fit__sky___handles_special_behaviour(masked_imaging_7x7):
+    masked_imaging_7x7 = copy.copy(masked_imaging_7x7)
+
+    masked_imaging_7x7.data -= 100.0
+
+    g0 = al.Galaxy(
+        redshift=0.5,
+        bulge=al.lp.Sersic(intensity=1.0),
+        disk=al.lp.Sersic(intensity=2.0),
+        mass_profile=al.mp.IsothermalSph(einstein_radius=1.0),
+    )
+
+    g1 = al.Galaxy(redshift=1.0, bulge=al.lp.Sersic(intensity=1.0))
+
+    tracer = al.Tracer(galaxies=[g0, g1])
+
+    fit = al.FitImaging(
+        dataset=masked_imaging_7x7,
+        tracer=tracer,
+        sky=al.lp_linear.Sky(),
+        settings_inversion=al.SettingsInversion(use_positive_only_solver=False),
+    )
+
+    assert fit.perform_inversion is True
+    assert fit.figure_of_merit == pytest.approx(-733125.35694344, 1.0e-4)
+
+    sky = fit.sky_linear_light_profiles_to_light_profiles
+
+    assert sky.light_profile_list[0].intensity == pytest.approx(-737.44550397, 1.0e-4)
+    assert sky.light_profile_list[1].intensity == pytest.approx(737.44550397, 1.0e-4)
+
+    fit = al.FitImaging(
+        dataset=masked_imaging_7x7,
+        tracer=tracer,
+        sky=al.lp_linear.Sky(),
+        settings_inversion=al.SettingsInversion(use_positive_only_solver=True),
+    )
+
+    assert fit.perform_inversion is True
+    assert fit.figure_of_merit == pytest.approx(-733125.3569434, 1.0e-4)
+
+    sky = fit.sky_linear_light_profiles_to_light_profiles
+
+    assert sky.light_profile_list[0].intensity == pytest.approx(0.0, 1.0e-4)
+    assert sky.light_profile_list[1].intensity == pytest.approx(1474.891048, 1.0e-4)
+
+    fit = al.FitImaging(
+        dataset=masked_imaging_7x7, tracer=tracer, sky=al.lp.Sky(intensity=-99.0)
+    )
+
+    assert fit.perform_inversion is False
+    assert fit.figure_of_merit == pytest.approx(-2862836.077500, 1.0e-4)
+
+
 def test__galaxy_model_image_dict(masked_imaging_7x7):
 
     # Normal Light Profiles Only