feat: make the minimum MGE Gaussian sigma configurable (sigma_min)

autogalaxy/analysis/model_util.py

@@ -17,6 +17,7 @@ def mge_model_from(
     ell_comps_uniform_width: float = 0.2,
     ell_comps_sigma : float = 0.3,
     use_spherical: bool = False,
+    sigma_min: float = 1e-4,
 ) -> af.Collection:
     """
     Construct a Multi-Gaussian Expansion (MGE) for the lens or source galaxy light.
@@ -81,6 +82,11 @@ def mge_model_from(
     use_spherical
         If True, use ``GaussianSph`` (no ell_comps). If False (default), use
         ``Gaussian`` with ellipticity.
+    sigma_min
+        The smallest Gaussian width (`sigma`) in arcseconds, which sets the lower end
+        of the log-spaced sigma values. Defaults to ``1e-4``. Increase it (e.g. to a
+        tenth of the pixel scale) to stop the basis wasting components on scales the
+        data cannot resolve.
 
     Returns
     -------
@@ -97,8 +103,22 @@ def mge_model_from(
     from autogalaxy.profiles.light.linear import Gaussian, GaussianSph
     from autogalaxy.profiles.basis import Basis
 
-    # The sigma values of the Gaussians will be fixed to values spanning 0.01 to the mask radius.
-    log10_sigma_list = np.linspace(-4, np.log10(mask_radius), total_gaussians)
+    if sigma_min <= 0.0:
+        raise ValueError(
+            f"mge_model_from requires sigma_min > 0.0, got {sigma_min}."
+        )
+
+    if sigma_min > mask_radius:
+        raise ValueError(
+            f"mge_model_from requires sigma_min <= mask_radius, got sigma_min="
+            f"{sigma_min} and mask_radius={mask_radius}."
+        )
+
+    # The sigma values of the Gaussians are fixed to log-spaced values spanning
+    # `sigma_min` (default 0.0001") to the mask radius.
+    log10_sigma_list = np.linspace(
+        np.log10(sigma_min), np.log10(mask_radius), total_gaussians
+    )
 
     if use_spherical:
         model_cls = GaussianSph
@@ -179,13 +199,14 @@ def mge_point_model_from(
     pixel_scales: float,
     total_gaussians: int = 10,
     centre: Tuple[float, float] = (0.0, 0.0),
+    sigma_min: float = 0.01,
 ) -> af.Model:
     """
     Construct a Multi-Gaussian Expansion (MGE) model for a compact or unresolved
     point-like component (e.g. a nuclear starburst, AGN, or unresolved bulge).
 
     The model is composed of ``total_gaussians`` linear Gaussians whose sigma values
-    are logarithmically spaced between 0.01 arcseconds and twice the pixel scale.
+    are logarithmically spaced between ``sigma_min`` and twice the pixel scale.
     All Gaussians share the same centre and ellipticity components, keeping the
     parameter count low while capturing a realistic PSF-convolved point source.
 
@@ -200,6 +221,11 @@ def mge_point_model_from(
     centre
         (y, x) centre of the point source in arc-seconds.  A ±0.1 arcsecond uniform
         prior is placed on each coordinate.
+    sigma_min
+        The smallest Gaussian width (`sigma`) in arcseconds, which sets the lower end
+        of the log-spaced sigma values. Defaults to ``0.01``. Increase it (e.g. to a
+        tenth of the pixel scale) to stop the basis wasting components on scales the
+        data cannot resolve.
 
     Returns
     -------
@@ -220,14 +246,19 @@ def mge_point_model_from(
             f"mge_point_model_from requires pixel_scales > 0, got {pixel_scales}."
         )
 
-    # Sigma values are logarithmically spaced between 0.01 arcsec (10**-2)
-    # and twice the pixel scale, with a floor to avoid taking log10 of
-    # very small or non-positive values.
-    min_log10_sigma = -2.0  # corresponds to 0.01 arcsec
-    max_sigma = max(2.0 * pixel_scales, 10**min_log10_sigma)
-    max_log10_sigma = np.log10(max_sigma)
+    if sigma_min <= 0.0:
+        raise ValueError(
+            f"mge_point_model_from requires sigma_min > 0.0, got {sigma_min}."
+        )
 
-    log10_sigma_list = np.linspace(min_log10_sigma, max_log10_sigma, total_gaussians)
+    # Sigma values are logarithmically spaced between `sigma_min` (default 0.01")
+    # and twice the pixel scale, with a floor to keep the upper end of the list
+    # at or above `sigma_min` when the pixel scale is very small.
+    max_sigma = max(2.0 * pixel_scales, sigma_min)
+
+    log10_sigma_list = np.linspace(
+        np.log10(sigma_min), np.log10(max_sigma), total_gaussians
+    )
     centre_0 = af.UniformPrior(lower_limit=centre[0] - 0.1, upper_limit=centre[0] + 0.1)
     centre_1 = af.UniformPrior(lower_limit=centre[1] - 0.1, upper_limit=centre[1] + 0.1)
 

test_autogalaxy/analysis/test_model_util.py

@@ -105,6 +105,72 @@ def test__mge_model_from__total_gaussians_per_basis():
     assert len(instance.profile_list) == 10
 
 
+def test__mge_model_from__sigma_min_default_spans_1e4_to_mask_radius():
+    model = ag.model_util.mge_model_from(mask_radius=3.0, total_gaussians=5)
+
+    instance = model.instance_from_prior_medians()
+    sigma_list = [profile.sigma for profile in instance.profile_list]
+
+    assert sigma_list[0] == pytest.approx(1e-4, 1.0e-8)
+    assert sigma_list[-1] == pytest.approx(3.0, 1.0e-8)
+
+
+def test__mge_model_from__sigma_min_input_sets_smallest_gaussian():
+    model = ag.model_util.mge_model_from(
+        mask_radius=3.0, total_gaussians=5, sigma_min=0.01
+    )
+
+    instance = model.instance_from_prior_medians()
+    sigma_list = [profile.sigma for profile in instance.profile_list]
+
+    assert sigma_list[0] == pytest.approx(0.01, 1.0e-8)
+    assert sigma_list[-1] == pytest.approx(3.0, 1.0e-8)
+    assert sigma_list == pytest.approx(
+        list(10 ** np.linspace(np.log10(0.01), np.log10(3.0), 5)), 1.0e-8
+    )
+
+
+def test__mge_model_from__sigma_min_invalid_raises():
+    with pytest.raises(ValueError):
+        ag.model_util.mge_model_from(
+            mask_radius=3.0, total_gaussians=5, sigma_min=0.0
+        )
+
+    with pytest.raises(ValueError):
+        ag.model_util.mge_model_from(
+            mask_radius=3.0, total_gaussians=5, sigma_min=4.0
+        )
+
+
+def test__mge_model_from__default_sigma_list_is_bitwise_unchanged():
+    """
+    The default `sigma_min=1e-4` must reproduce the hardcoded `np.linspace(-4, ...)`
+    ladder that predates the `sigma_min` argument EXACTLY, not approximately.
+
+    Every fixed `sigma` feeds the PyAutoFit identifier of a run, so drift here gives
+    existing fits a new `unique_id`, orphaning their output directories and silently
+    restarting them from scratch. The identifier quantizes floats at
+    `RESOLUTION = 1e-8` (see `autofit.mapper.identifier`), so it does not in fact move
+    for drift below that -- but exact equality is the stronger guarantee and costs
+    nothing, catching drift ~8 orders of magnitude earlier than the identifier does.
+
+    `pytest.approx(rel=1e-8)` is deliberately NOT used: it only fails once the ladder
+    has moved by a relative ~1e-7, which is already past the point where the
+    identifier changes.
+    """
+    for mask_radius, total_gaussians in [(3.0, 20), (3.5, 30), (7.5, 10), (1.0, 5)]:
+        model = ag.model_util.mge_model_from(
+            mask_radius=mask_radius, total_gaussians=total_gaussians
+        )
+
+        instance = model.instance_from_prior_medians()
+        sigma_list = [profile.sigma for profile in instance.profile_list]
+
+        assert sigma_list == list(
+            10 ** np.linspace(-4, np.log10(mask_radius), total_gaussians)
+        )
+
+
 def test__mge_point_model_from__returns_basis_model_with_correct_gaussians():
     """
     mge_point_model_from should return an af.Model wrapping a Basis whose
@@ -136,6 +202,51 @@ def test__mge_point_model_from__sigma_values_span_correct_range():
     assert gaussian_list[-1].sigma == pytest.approx(pixel_scales * 2.0, rel=1.0e-4)
 
 
+def test__mge_point_model_from__default_sigma_list_is_bitwise_unchanged():
+    """
+    As for `mge_model_from`, the default `sigma_min=0.01` must reproduce the
+    hardcoded `min_log10_sigma = -2.0` ladder that predates the argument EXACTLY,
+    so the identifier of an existing point-source fit does not change.
+    """
+    for pixel_scales, total_gaussians in [(0.1, 10), (0.05, 5), (0.2, 3), (0.001, 4)]:
+        model = ag.model_util.mge_point_model_from(
+            pixel_scales=pixel_scales, total_gaussians=total_gaussians
+        )
+
+        sigma_list = [gaussian.sigma for gaussian in model.profile_list]
+
+        max_sigma = max(2.0 * pixel_scales, 10**-2.0)
+
+        assert sigma_list == list(
+            10 ** np.linspace(-2.0, np.log10(max_sigma), total_gaussians)
+        )
+
+
+def test__mge_point_model_from__sigma_min_input_sets_smallest_gaussian():
+    total_gaussians = 5
+
+    model = ag.model_util.mge_point_model_from(
+        pixel_scales=0.1, total_gaussians=total_gaussians, sigma_min=0.01 / 10.0
+    )
+
+    sigma_list = [gaussian.sigma for gaussian in model.profile_list]
+
+    assert sigma_list[0] == pytest.approx(0.001, 1.0e-8)
+    assert sigma_list[-1] == pytest.approx(0.2, 1.0e-8)
+    assert sigma_list == pytest.approx(
+        list(10 ** np.linspace(np.log10(0.001), np.log10(0.2), total_gaussians)),
+        1.0e-8,
+    )
+
+
+def test__mge_point_model_from__sigma_min_invalid_raises():
+    with pytest.raises(ValueError):
+        ag.model_util.mge_point_model_from(pixel_scales=0.1, sigma_min=0.0)
+
+    with pytest.raises(ValueError):
+        ag.model_util.mge_point_model_from(pixel_scales=0.1, sigma_min=-1.0)
+
+
 def test__mge_point_model_from__shared_centre_and_ell_comps():
     """
     All Gaussians must share exactly the same centre prior objects and ell_comps