Skip to content
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
149 changes: 85 additions & 64 deletions autolens/analysis/model_util.py
Original file line number Diff line number Diff line change
Expand Up @@ -13,88 +13,109 @@
- ``hilbert_pixels_from_pixel_scale`` — estimate Hilbert image-mesh pixel count.

PyAutoLens-specific:
- ``simulator_start_here_model_from`` — builds a default imaging simulator model with
optional lens light and point-source components for start_here scripts.
- ``random_galaxies_for_simulation_from`` — sample concrete (lens, source) ``Galaxy``
instances for synthetic-data generation in ``start_here`` scripts.
"""
import autofit as af
from typing import Optional, Tuple

import numpy as np

import autolens as al

from autogalaxy.analysis.model_util import mge_model_from
from autogalaxy.analysis.model_util import mge_point_model_from
from autogalaxy.analysis.model_util import hilbert_pixels_from_pixel_scale

def simulator_start_here_model_from(
include_lens_light: bool = True, use_point_source: bool = False
):

if include_lens_light:
bulge = af.Model(al.lp_snr.Sersic)
SIMULATOR_RANDOM_LENS_SUMMARY = (
"Each simulated strong lens draws fresh truths from: "
"lens bulge SNR in [20, 60] (when included), "
"lens mass einstein_radius in [0.2, 1.8] with normal-clipped ellipticity, "
"external shear ~ Normal(0, 0.05), "
"source bulge SNR in [10, 30] / point-source flux in [0.0, 2.0] (mode dependent)."
)


def _clipped_ell_comp(rng: np.random.Generator) -> float:
return float(np.clip(rng.normal(0.0, 0.2), -1.0, 1.0))


def random_galaxies_for_simulation_from(
include_lens_light: bool = True,
use_point_source: bool = False,
rng: Optional[np.random.Generator] = None,
) -> Tuple["al.Galaxy", "al.Galaxy"]:
"""
Sample a ``(lens_galaxy, source_galaxy)`` pair for synthetic strong-lens
data generation.

Each parameter is drawn directly from a numpy ``Generator`` and used to
construct concrete profile instances — no ``af.Model`` priors are involved.
SNR-normalised Sersic profiles (``lp_snr.Sersic``) are used for diffuse
light components so that simulator output lands at a controlled target
SNR; the SNR appears as a profile attribute on the *instance*, never as a
fitting parameter.

Do **not** use the returned galaxies as fitting models. They are
instances, suitable for ``Tracer`` / ``simulator.via_tracer_from``.

Parameters
----------
include_lens_light
If True (default), give the lens galaxy an ``lp_snr.Sersic`` bulge.
If False, the lens is mass-only.
use_point_source
If True, source is a ``PointFlux`` with random centre and flux. If
False (default), source is an ``lp_snr.Sersic``.
rng
Optional ``numpy.random.Generator``. If ``None`` a fresh
``default_rng()`` is created on each call.

Returns
-------
(Galaxy, Galaxy)
``(lens_galaxy, source_galaxy)`` at redshifts 0.5 and 1.0 respectively.
"""
rng = rng if rng is not None else np.random.default_rng()

bulge.centre = (0.0, 0.0)
bulge.ell_comps.ell_comps_0 = af.TruncatedGaussianPrior(
mean=0.0, sigma=0.2, lower_limit=-1.0, upper_limit=1.0
)
bulge.ell_comps.ell_comps_1 = af.TruncatedGaussianPrior(
mean=0.0, sigma=0.2, lower_limit=-1.0, upper_limit=1.0
)
bulge.signal_to_noise_ratio = af.UniformPrior(
lower_limit=20.0, upper_limit=60.0
)
bulge.effective_radius = af.UniformPrior(lower_limit=1.0, upper_limit=5.0)
bulge.sersic_index = af.TruncatedGaussianPrior(
mean=4.0, sigma=0.5, lower_limit=0.8, upper_limit=5.0
if include_lens_light:
lens_bulge = al.lp_snr.Sersic(
centre=(0.0, 0.0),
ell_comps=(_clipped_ell_comp(rng), _clipped_ell_comp(rng)),
effective_radius=float(rng.uniform(1.0, 5.0)),
sersic_index=float(rng.uniform(3.5, 4.5)),
signal_to_noise_ratio=float(rng.uniform(20.0, 60.0)),
)
else:
bulge = None

mass = af.Model(al.mp.Isothermal)
lens_bulge = None

mass.centre = (0.0, 0.0)
mass.ell_comps.ell_comps_0 = af.TruncatedGaussianPrior(
mean=0.0, sigma=0.2, lower_limit=-1.0, upper_limit=1.0
)
mass.ell_comps.ell_comps_1 = af.TruncatedGaussianPrior(
mean=0.0, sigma=0.2, lower_limit=-1.0, upper_limit=1.0
mass = al.mp.Isothermal(
centre=(0.0, 0.0),
ell_comps=(_clipped_ell_comp(rng), _clipped_ell_comp(rng)),
einstein_radius=float(rng.uniform(0.2, 1.8)),
)
mass.einstein_radius = af.UniformPrior(lower_limit=0.2, upper_limit=1.8)

shear = af.Model(al.mp.ExternalShear)

shear.gamma_1 = af.GaussianPrior(mean=0.0, sigma=0.05)
shear.gamma_2 = af.GaussianPrior(mean=0.0, sigma=0.05)
shear = al.mp.ExternalShear(
gamma_1=float(rng.normal(0.0, 0.05)),
gamma_2=float(rng.normal(0.0, 0.05)),
)

lens = af.Model(al.Galaxy, redshift=0.5, bulge=bulge, mass=mass, shear=shear)
lens = al.Galaxy(redshift=0.5, bulge=lens_bulge, mass=mass, shear=shear)

if use_point_source:

point_0 = af.Model(al.ps.PointFlux)

point_0.centre.centre_0 = af.GaussianPrior(mean=0.0, sigma=0.3)
point_0.centre.centre_1 = af.GaussianPrior(mean=0.0, sigma=0.3)
point_0.flux = af.UniformPrior(lower_limit=0.0, upper_limit=2.0)

source = af.Model(al.Galaxy, redshift=1.0, point_0=point_0)

else:

bulge = af.Model(al.lp_snr.Sersic)

bulge.centre_0 = af.GaussianPrior(mean=0.0, sigma=0.3)
bulge.centre_1 = af.GaussianPrior(mean=0.0, sigma=0.3)
bulge.ell_comps.ell_comps_0 = af.TruncatedGaussianPrior(
mean=0.0, sigma=0.2, lower_limit=-1.0, upper_limit=1.0
)
bulge.ell_comps.ell_comps_1 = af.TruncatedGaussianPrior(
mean=0.0, sigma=0.2, lower_limit=-1.0, upper_limit=1.0
)
bulge.signal_to_noise_ratio = af.UniformPrior(
lower_limit=10.0, upper_limit=30.0
point_0 = al.ps.PointFlux(
centre=(float(rng.normal(0.0, 0.3)), float(rng.normal(0.0, 0.3))),
flux=float(rng.uniform(0.0, 2.0)),
)
bulge.effective_radius = af.UniformPrior(lower_limit=0.01, upper_limit=3.0)
bulge.sersic_index = af.TruncatedGaussianPrior(
mean=2.0, sigma=0.5, lower_limit=0.8, upper_limit=5.0
source = al.Galaxy(redshift=1.0, point_0=point_0)
else:
source_bulge = al.lp_snr.Sersic(
centre=(float(rng.normal(0.0, 0.3)), float(rng.normal(0.0, 0.3))),
ell_comps=(_clipped_ell_comp(rng), _clipped_ell_comp(rng)),
effective_radius=float(rng.uniform(0.01, 3.0)),
sersic_index=float(rng.uniform(1.5, 2.5)),
signal_to_noise_ratio=float(rng.uniform(10.0, 30.0)),
)

source = af.Model(al.Galaxy, redshift=1.0, bulge=bulge)
source = al.Galaxy(redshift=1.0, bulge=source_bulge)

return lens, source
Loading