Skip to content

Latest commit

 

History

History
946 lines (596 loc) · 34.6 KB

File metadata and controls

946 lines (596 loc) · 34.6 KB

✏️ This page is auto-generated from scripts/chapter_1_introduction/tutorial_1_grids_and_galaxies.py — do not edit it directly. It shows the example fully executed, with its real output images. Run it yourself via the Python script or the Jupyter notebook.

HowToLens: Introduction

A strong gravitational lens is a system where two (or more) galaxies align perfectly down our line of sight from Earth such that the foreground galaxy's mass curves space-time in on itself, such that the light of a background source galaxy is deflected and magnified. This means we can see the background source galaxy multiple times, as multiple arcs or rings, because multiple paths through the foreground galaxy's mass are taken by the source's light.

Here is a schematic of a strong gravitational lens:

Schematic of Gravitational Lensing Credit: F. Courbin, S. G. Djorgovski, G. Meylan, et al., Caltech / EPFL / WMKO https://www.astro.caltech.edu/~george/qsolens/

Today, Astronomers use computer software, statistical algorithms and image processing techniques analyse the light of these strong gravitational lenses. They are used for a wide range of scientific studies, including studying dark matter in the foreground lens galaxies, the morphology and structure of the background source galaxies, and even the expansion of the Universe itself.

The HowToLens series of tutorials will teach you how to perform this analysis yourself, using the open-source software package PyAutoLens. By the end of the HowToLens series, you'll be able to take an image of strong lens and study it using the same techniques that professional astronomers use today.

Tutorial 1: Grids And Galaxies

In this tutorial, we will introduce the first fundamental concepts and quantities used to study strong lenses. These concepts will enable us to create images of galaxies and analyze how their light is distributed across space. Additionally, we will explore how adjusting various properties of galaxies can alter their appearance. For instance, we can change the size of a galaxy, rotate it, or modify its brightness.

To create these images, we first need to define 2D grids of ((y, x)) coordinates. We will shift and rotate these grids to manipulate the appearance of the galaxy in the generated images. The grid will serve as the input for light profiles, which are analytic functions that describe the distribution of a galaxy's light. By evaluating these light profiles on the grid, we can effectively generate images that represent the structure and characteristics of galaxies.

This tutorial won't yet perform any lensing calculations, which are introduced in the next tutorial.

Here is an overview of what we'll cover in this tutorial:

  • Grids: We'll create a uniform grid of $(y,x)$ coordinates and show how it can be used to measure the light of a galaxy.
  • Geometry: We'll show how to shift and rotate a grid, and convert it to elliptical coordinates.
  • Light Profiles: We'll introduce light profiles, analytic functions that describe how a galaxy's light is distributed.
  • Galaxies: We'll create galaxies containing light profiles and show how to compute their image.
  • Units: We'll show how to convert the units of a galaxy's image to physical units like kiloparsecs.

The imports below are required to run the HowToLens tutorials in a Jupiter notebook. They also import the autolens package and the autolens.plot module which are used throughout the tutorials.

Contents

  • Grids: A Grid2D is a set of two-dimensional $(y,x)$ coordinates that represent points in space where we.
  • Geometry: The above grid is centered on the origin (0.0", 0.0").
  • Light Profiles: Galaxies are collections of stars, gas, dust, and other astronomical objects that emit light.
  • One Dimension Projection: We often want to calculative 1D quantities of a light profile, for example to plot how its light.
  • Galaxies: Galaxies are collections of light profiles that represent a galaxy's luminous emission.
  • Units: By assuming a redshift for a galaxy we can convert its quantities from arcseconds to kiloparsecs.
  • Wrap Up: Summary of the script and next steps.
  • Advanced Topics: Optional reference material for expert users, not needed on a first read.
  • Other Unit Conversion: Other unit conversions a galaxy's redshift enables (e.g. AB magnitudes).
  • Over Sampling: The adaptive over-sampling scheme used to evaluate light profiles accurately.
from autolens import jax_wrapper  # Sets JAX environment before other imports

from autolens import setup_notebook; setup_notebook()

import matplotlib.pyplot as plt
import numpy as np

import autoarray as aa
import autolens as al
import autolens.plot as aplt
.../PyAutoNerves/autonerves/workspace.py:206: UserWarning: Cannot verify the workspace at HowToLens/scripts/chapter_1_introduction is compatible with the installed library version (2026.7.23.1): no `version.minimum_library_version` or `version.workspace_version` key in config/general.yaml and no version.txt at the workspace root.

If you cloned the workspace from `main` rather than a release tag, set `version.workspace_version_check: False` in config/general.yaml to silence this warning. The `main` branch updates more frequently than library releases, so version mismatches are expected and not actionable for `main`-branch users.

You can also set the environment variable PYAUTO_SKIP_WORKSPACE_VERSION_CHECK=1 to disable temporarily.
  warnings.warn(_missing_version_warning(root, library_version))
.../PyAutoNerves/autonerves/workspace.py:206: UserWarning: Cannot verify the workspace at HowToLens/scripts/chapter_1_introduction is compatible with the installed library version (2026.7.23.1): no `version.minimum_library_version` or `version.workspace_version` key in config/general.yaml and no version.txt at the workspace root.

If you cloned the workspace from `main` rather than a release tag, set `version.workspace_version_check: False` in config/general.yaml to silence this warning. The `main` branch updates more frequently than library releases, so version mismatches are expected and not actionable for `main`-branch users.

You can also set the environment variable PYAUTO_SKIP_WORKSPACE_VERSION_CHECK=1 to disable temporarily.
  warnings.warn(_missing_version_warning(root, library_version))
Working Directory has been set to `HowToLens`
.../PyAutoNerves/autonerves/workspace.py:206: UserWarning: Cannot verify the workspace at HowToLens/scripts/chapter_1_introduction is compatible with the installed library version (2026.7.23.1): no `version.minimum_library_version` or `version.workspace_version` key in config/general.yaml and no version.txt at the workspace root.

If you cloned the workspace from `main` rather than a release tag, set `version.workspace_version_check: False` in config/general.yaml to silence this warning. The `main` branch updates more frequently than library releases, so version mismatches are expected and not actionable for `main`-branch users.

You can also set the environment variable PYAUTO_SKIP_WORKSPACE_VERSION_CHECK=1 to disable temporarily.
  warnings.warn(_missing_version_warning(root, library_version))

Grids

A Grid2D is a set of two-dimensional $(y,x)$ coordinates that represent points in space where we evaluate the light emitted by a galaxy.

Each coordinate on the grid is referred to as a 'pixel'. This is because we use the grid to measure the brightness of a galaxy at each of these coordinates, allowing us to create an image of the galaxy.

Grids are defined in units of 'arc-seconds' ("). An arc-second is a unit of angular measurement used by astronomers to describe the apparent size of objects in the sky.

The pixel_scales parameter sets how many arc-seconds each pixel represents. For example, if pixel_scales=0.1, then each pixel covers 0.1" of the sky.

We create a uniform 2D grid of 101 x 101 pixels with a pixel scale of 0.1", corresponding to an area of 10.1" x 10.1", spanning from -5.05" to 5.05" in both the y and x directions.

grid = al.Grid2D.uniform(
    shape_native=(
        101,
        101,
    ),  # The dimensions of the grid, which here is 101 x 101 pixels.
    pixel_scales=0.1,  # The conversion factor between pixel units and arc-seconds.
)

We can visualize this grid as a uniform grid of dots, each representing a coordinate where the light is measured.

aplt.plot_grid(grid=grid, title="Uniform Grid of Coordinates")

png

Each coordinate in the grid corresponds to an arc-second position. Below, we print a few of these coordinates to see the values.

print("(y,x) pixel 0:")
print(grid.native[0, 0])  # The coordinate of the first pixel.
print("(y,x) pixel 1:")
print(grid.native[0, 1])  # The coordinate of the second pixel.
print("(y,x) pixel 2:")
print(grid.native[0, 2])  # The coordinate of the third pixel.
print("(y,x) pixel 100:")
print(grid.native[1, 0])  # The coordinate of the 100th pixel.
print("...")
(y,x) pixel 0:
[ 5. -5.]
(y,x) pixel 1:
[ 5.  -4.9]
(y,x) pixel 2:
[ 5.  -4.8]
(y,x) pixel 100:
[ 4.9 -5. ]
...

Grids have two internal representations, native and slim:

  • native: A 2D array with shape [total_y_pixels, total_x_pixels, 2], where the 2 corresponds to the (y,x) coordinates.
  • slim: A 1D array with shape [total_y_pixels * total_x_pixels, 2], where the coordinates are 'flattened' into a single list.

These formats are useful for different calculations and plotting. Here, we show the same coordinate using both formats.

print("(y,x) pixel 0 (accessed via native):")
print(grid.native[0, 0])
print("(y,x) pixel 0 (accessed via slim 1D):")
print(grid.slim[0])
(y,x) pixel 0 (accessed via native):
[ 5. -5.]
(y,x) pixel 0 (accessed via slim 1D):
[ 5. -5.]

We can also check the shapes of the Grid2D object in both native and slim formats. For this grid, the native shape is (101, 101, 2) and the slim shape is (10201, 2).

print(grid.native.shape)
print(grid.slim.shape)
(101, 101, 2)
(10201, 2)

For the HowToLens tutorials, you don't need to fully understand why grids have both native and slim representations. Just note that both are used for calculations and plotting.

Exercise: Try creating grids with different shapes and pixel scales using the al.Grid2D.uniform() function above. Observe how the grid coordinates change when you adjust shape_native and pixel_scales.

Geometry

The above grid is centered on the origin (0.0", 0.0"). Sometimes, we need to shift the grid to be centered on a specific point, like the center of a galaxy.

We can shift the grid to a new center, (y_c, x_c), by subtracting this center from each coordinate.

centre = (0.3, 0.5)  # Shifting the grid to be centered at y=0.3", x=0.5".

grid_shifted = grid
grid_shifted[:, 0] = grid_shifted[:, 0] - centre[0]  # Shift in y-direction.
grid_shifted[:, 1] = grid_shifted[:, 1] - centre[1]  # Shift in x-direction.

print("(y,x) pixel 0 After Shift:")
print(grid_shifted.native[0, 0])  # The coordinate of the first pixel after shifting.
(y,x) pixel 0 After Shift:
[ 4.7 -5.5]

The grid is now centered around (0.3", 0.5"). We can plot the shifted grid to see this change.

Exercise: Try shifting the grid to a different center, for example (0.0", 0.0") or (2.0", 3.0"). Observe how the center of the grid changes when you adjust the centre variable.

aplt.plot_grid(grid=grid_shifted, title="Grid Centered Around (0.3, 0.5)")

png

Next, we can rotate the grid by an angle phi (in degrees). The rotation is counter-clockwise from the positive x-axis.

To rotate the grid:

  1. Calculate the distance radius of each coordinate from the origin using $r = \sqrt{y^2 + x^2}$.
  2. Determine the angle theta counter clockwise from the positive x-axis using $\theta = \arctan(y / x)$.
  3. Adjust theta by the rotation angle and convert back to Cartesian coordinates via $y = r \sin(\theta)$ and $x = r \cos(\theta)$.
angle_degrees = 60.0

y = grid_shifted[:, 0]
x = grid_shifted[:, 1]

radius = np.sqrt(y**2 + x**2)
theta = np.arctan2(y, x) - np.radians(angle_degrees)

grid_rotated = grid_shifted
grid_rotated[:, 0] = radius * np.sin(theta)
grid_rotated[:, 1] = radius * np.cos(theta)

print("(y,x) pixel 0 After Rotation:")
print(grid_rotated.native[0, 0])  # The coordinate of the first pixel after rotation.
(y,x) pixel 0 After Rotation:
[7.11313972 1.3203194 ]

The grid has now been rotated 60 degrees counter-clockwise. We can plot it to see the change.

Exercise: Try rotating the grid by a different angle, for example 30 degrees or 90 degrees. Observe how the grid changes when you adjust the angle_degrees variable.

aplt.plot_grid(grid=grid_rotated, title="Grid Rotated 60 Degrees")

png

Next, we convert the rotated grid to elliptical coordinates using:

$\eta = \sqrt{(x_r)^2 + (y_r)^2/q^2}$

Where q is the axis-ratio of the ellipse and (x_r, y_r) are the rotated coordinates.

Elliptical coordinates are a system used to describe positions in relation to an ellipse rather than a circle. They are particularly useful in astronomy when dealing with objects like galaxies, which often have elliptical shapes due to their inclination or intrinsic shape.

Exercise: Try converting the grid to elliptical coordinates using a different axis-ratio, for example 0.3 or 0.8. What happens to the grid when you adjust the axis_ratio variable?

axis_ratio = 0.5
eta = np.sqrt((grid_rotated[:, 0]) ** 2 + (grid_rotated[:, 1]) ** 2 / axis_ratio**2)

Above, the angle $\phi$ (in degrees) was used to rotate the grid, and the axis-ratio $q$ was used to convert the grid to elliptical coordinates.

From now on, we'll describe ellipticity using "elliptical components" $\epsilon_{1}$ and $\epsilon_{2}$, calculated from $\phi$ and $q$:

$\epsilon_{1} = \frac{1 - q}{1 + q} \sin(2\phi)$
$\epsilon_{2} = \frac{1 - q}{1 + q} \cos(2\phi)$

We'll refer to these as ell_comps in the code for brevity.

Future tutorials will explain why $\epsilon_{1}$ and $\epsilon_{2}$ are preferred over $q$ and $\phi$.

Exercise: Try computing the elliptical components from the axis-ratio and angle above. What happens to the elliptical components when you adjust the axis_ratio and angle_degrees variables?

fac = (1 - axis_ratio) / (1 + axis_ratio)
epsilon_y = fac * np.sin(2 * np.radians(angle_degrees))
epsilon_x = fac * np.cos(2 * np.radians(angle_degrees))

ell_comps = (epsilon_y, epsilon_x)

print("Elliptical Components:")
print(ell_comps)
Elliptical Components:
(np.float64(0.28867513459481287), np.float64(-0.16666666666666657))

Light Profiles

Galaxies are collections of stars, gas, dust, and other astronomical objects that emit light. Astronomers study this light to understand various properties of galaxies.

To model the light of a galaxy, we use light profiles, which are mathematical functions that describe how a galaxy's light is distributed across space. By applying these light profiles to 2D grids of $(y, x)$ coordinates, we can create images that represent a galaxy's luminous emission.

A commonly used light profile is the Sersic profile, which is widely adopted in astronomy for representing galaxy light. The Sersic profile is defined by the equation:

$I_{\rm Ser} (\eta_{\rm l}) = I \exp \left{ -k \left[ \left( \frac{\eta}{R} \right)^{\frac{1}{n}} - 1 \right] \right}$

In this equation:

  • $\eta$ represents the elliptical coordinates of the profile in arc-seconds (refer to earlier sections for elliptical coordinates).
  • $I$ is the intensity normalization of the profile, given in arbitrary units, which controls the overall brightness of the Sersic profile.
  • $R$ is the effective radius in arc-seconds, which determines the size of the profile.
  • $n$ is the Sersic index, which defines how 'steep' the profile is, influencing the concentration of light.
  • $k$ is a constant that ensures half the light of the profile lies within the radius $R$, where $k = 2n - \frac{1}{3}$.

We can evaluate this function using values for $(\eta, I, R, n)$ to calculate the intensity of the profile at a particular elliptical coordinate.

elliptical_coordinate = (
    0.5  # The elliptical coordinate where we compute the intensity, in arc-seconds.
)
intensity = 1.0  # Intensity normalization of the profile in arbitrary units.
effective_radius = 2.0  # Effective radius of the profile in arc-seconds.
sersic_index = 1.0  # Sersic index of the profile.
k = 2 * sersic_index - (
    1.0 / 3.0
)  # Calculating the constant k, note that this is an approximation.

# Calculate the intensity of the Sersic profile at a specific elliptical coordinate.
sersic_value = np.exp(
    -k * ((elliptical_coordinate / effective_radius) ** (1.0 / sersic_index) - 1.0)
)

print("Intensity of Sersic Light Profile at Elliptical Coordinate 0.5:")
print(sersic_value)
Intensity of Sersic Light Profile at Elliptical Coordinate 0.5:
3.4903429574618414

The calculation above gives the intensity of the Sersic profile at an elliptical coordinate of 0.5.

To create a complete image of the Sersic profile, we can evaluate the intensity at every point in our grid of elliptical coordinates.

sersic_image = np.exp(-k * ((eta / effective_radius) ** (1.0 / sersic_index) - 1.0))

When we plot the resulting image, we can see how the properties of the grid affect its appearance:

  • The peak intensity is at the position (0.3", 0.5"), where we shifted the grid.
  • The image is elongated along a 60° counter-clockwise angle, corresponding to the rotation of the grid.
  • The image has an elliptical shape, consistent with the axis ratio of 0.5.

This demonstrates how the geometry of the grid directly influences the appearance of the light profile.

Exercise: Try changing the values of centre, ell_comps, effective_radius, and sersic_index above. Observe how these adjustments change the Sersic profile image.

aplt.plot_array(
    array=aa.Array2D(values=sersic_image, mask=grid.mask), title="Sersic Image"
)

png

Instead of manually handling these transformations, we can use LightProfile objects from the light_profile module (lp) for faster and more efficient calculations.

Below, we define a Sersic light profile using the Sersic object. We can print the profile to display its parameters.

sersic_light_profile = al.lp.Sersic(
    centre=(0.0, 0.0),
    ell_comps=(0.0, 0.1),
    intensity=1.0,
    effective_radius=2.0,
    sersic_index=1.0,
)

print(sersic_light_profile)
Sersic
centre: (0.0, 0.0)
ell_comps: (0.0, 0.1)
intensity: 1.0
effective_radius: 2.0
sersic_index: 1.0

With this Sersic light profile, we can create an image by passing a grid to its image_2d_from method.

The calculation will internally handle all the coordinate transformations and intensity evaluations we performed manually earlier, making it much simpler.

The Sersic profile we created just above is different from the one we used to manually compute the image, so the image will look different. However, the process is the same.

image = sersic_light_profile.image_2d_from(grid=grid)

aplt.plot_array(array=image, title="Sersic Image via Light Profile")

png

The image is returned as an Array2D object. Similar to a Grid2D, it has two forms:

  • native: A 2D array with shape [total_y_image_pixels, total_x_image_pixels].
  • slim: A 1D array that flattens this data into shape [total_y_image_pixels * total_x_image_pixels].

The native form is often used for visualizations, while the slim form can be useful for certain calculations.

print("Intensity of pixel 0:")
print(image.native[0, 0])
print("Intensity of pixel 1:")
print(image.slim[1])
Intensity of pixel 0:
0.013366491142097661
Intensity of pixel 1:
0.014020623587576735

To visualize the light profile's image, we use aplt.plot_array.

We provide it with the light profile and the grid, which are used to create and plot the image.

aplt.plot_array(
    array=sersic_light_profile.image_2d_from(grid=grid), title="Image via Light Profile"
)

png

One Dimension Projection

We often want to calculative 1D quantities of a light profile, for example to plot how its light changes as a function of radius.

To do this, we must still input a 2D grid into the image_2d_from method, therefore we create a project 2D radial grid as follows which has shape [Number_of_1d_coordinates, 2] and where all [:,0] entries are the same.

A simple example of such a grid is as follows with 4 1D coordinates is:

grid_2d_projected = al.Grid2DIrregular(
    [
        [1.000000e-06, 1.000000e-06],
        [1.000000e-06, 1.000001e00],
        [1.000000e-06, 2.000001e00],
        [1.000000e-06, 3.000001e00],
    ]
)

As in this example, we often already have a 2D grid we are using to calculate images of a ligth profile and it would be convenient to simply create grid_2d_projected from that.

For example, we may want the project grid which traces it major axis in uniform radial steps.

This is easily computed using the grid_2d_radial_project_from function and passing the centre and angle of a light profile we can make it align with the light profile itself.

Note how in this example the two galaxy bulges are not rotationally aligned but we aligned the projected grid with the first galaxy. The centres are aligned, but if they were not that would cause similar issues.

grid_2d_projected = grid.grid_2d_radial_projected_from(
    centre=sersic_light_profile.centre, angle=sersic_light_profile.angle()
)

image_1d = sersic_light_profile.image_2d_from(grid=grid_2d_projected)

We can now plot the 1D radial profile of the light profile. This profile shows how the intensity of the light changes as a function of distance from the profile's center. This is a more informative way to visualize the light p rofile's distribution.

When we plot 1D quantities, we do not use built-in plotting functions as in 2D, but instead use standard matplotlib functionality.

The reason is partly that 1D plotting is simple, but also because 1D plots have many different decisions about what is plotted and how they are computed, meaning its better to give the user full control.

Exercise: Try plotting the 1D radial profile of Sersic profiles with different effective radii and Sersic indices. Does the 1D representation show more clearly how the light distribution changes with these parameters?

plt.plot(grid_2d_projected[:, 1], image_1d)
plt.xlabel("Radius (arcseconds)")
plt.ylabel("Luminosity")
plt.show()
plt.close()

png

Since galaxy light distributions often cover a wide range of values, they are typically better visualized on a log10 scale. This approach helps highlight details in the faint outskirts of a light profile.

The plot_array/subplot_\* object has a use_log10 option that applies this transformation automatically. Below, you can see that the image plotted in log10 space reveals more details.

aplt.plot_array(
    array=sersic_light_profile.image_2d_from(grid=grid),
    title="Sersic Image",
    use_log10=True,
)

png

Galaxies

Now, let's introduce Galaxy objects, which are a key component in PyAutoLens.

A light profile represents a single feature of a galaxy, such as its bulge or disk. To model a complete galaxy, we combine multiple light profiles into a Galaxy object. This allows us to create images that include different components of a galaxy.

In addition to light profiles, a Galaxy has a redshift, which indicates how far away it is from Earth. The redshift is essential for performing unit conversions using cosmological calculations, such as converting arc-seconds into kiloparsecs. (A kiloparsec is a distance unit in astronomy, equal to about 3.26 million light-years.)

Redshifts are especially important in strong lensing, where the foreground lens galaxy and background source galaxy lie at two different redshifts. We are not yet performing any lensing calculations in this tutorial, so for now we simply use a single galaxy to build up intuition for the Galaxy object.

Let's start by creating a galaxy with two Sersic light profiles, which we will consider to represent a bulge and disk component of the galaxy, the two most important structures seen in galaxies.

bulge = al.lp.Sersic(
    centre=(0.0, 0.0),
    ell_comps=(0.0, 0.111111),
    intensity=1.0,
    effective_radius=1.0,
    sersic_index=2.5,
)

disk = al.lp.Sersic(
    centre=(0.0, 0.0),
    ell_comps=(0.0, 0.3),
    intensity=0.3,
    effective_radius=3.0,
    sersic_index=1.0,
)

galaxy = al.Galaxy(redshift=0.5, bulge=bulge, disk=disk)

print(galaxy)
Redshift: 0.5
Light Profiles:
Sersic
centre: (0.0, 0.0)
ell_comps: (0.0, 0.111111)
intensity: 1.0
effective_radius: 1.0
sersic_index: 2.5
Sersic
centre: (0.0, 0.0)
ell_comps: (0.0, 0.3)
intensity: 0.3
effective_radius: 3.0
sersic_index: 1.0

We can pass a 2D grid to a light profile to compute its image using the image_2d_from method.

The same approach works for a Galaxy object:

image = galaxy.image_2d_from(grid=grid)

print("Intensity of `Grid2D` pixel 0:")
print(image.native[0, 0])
print("Intensity of `Grid2D` pixel 1:")
print(image.native[0, 1])
print("Intensity of `Grid2D` pixel 2:")
print(image.native[0, 2])
print("...")
Intensity of `Grid2D` pixel 0:
0.02489491716484805
Intensity of `Grid2D` pixel 1:
0.025428546280541576
Intensity of `Grid2D` pixel 2:
0.02596640780160061
...

We can plot the galaxy's image, just like how we did for a light profile.

aplt.plot_array(array=galaxy.image_2d_from(grid=grid), title="Galaxy Bulge+Disk Image")

png

The bulge dominates the center of the image, and is pretty much the only luminous emission we can see on a linear scale. The disk's emission is present, but it is much fainter and spread over a larger area.

We can confirm this using the subplot_galaxy_light_profiles method, which plots each individual light profile separately.

aplt.subplot_galaxy_light_profiles(galaxy=galaxy, grid=grid)

png

Because galaxy light distributions often follow a log10 pattern, plotting in log10 space helps reveal details in the outskirts of the light profile, in this case the emission of the disk.

This is especially helpful to separate the bulge and disk profiles, which have different intensities and sizes.

aplt.plot_array(
    array=galaxy.image_2d_from(grid=grid),
    title="Galaxy Bulge+Disk Image",
    use_log10=True,
)

png

Using the tools above, we can visualize each light profile's contribution in 1D.

1D plots show the intensity of the light profile as a function of distance from the profile's center. The bulge and disk profiles in this example share the same centre, meaning that plotting them together on the same 1D plot shows how they vary relative to one another.

If the centre of the profiles were different, when you make the 1D plot you would need to decide whether to plot the profiles offset from one another or plot them both from zero.

grid_2d_projected = grid.grid_2d_radial_projected_from(
    centre=galaxy.bulge.centre, angle=galaxy.bulge.angle()
)
bulge_image_1d = galaxy.bulge.image_2d_from(grid=grid_2d_projected)

grid_2d_projected = grid.grid_2d_radial_projected_from(
    centre=galaxy.disk.centre, angle=galaxy.disk.angle()
)
disk_image_1d = galaxy.disk.image_2d_from(grid=grid_2d_projected)

plt.plot(grid_2d_projected[:, 1], bulge_image_1d, label="Bulge")
plt.plot(grid_2d_projected[:, 1], disk_image_1d, label="Disk")
plt.xlabel("Radius (arcseconds)")
plt.ylabel("Luminosity")
plt.legend()
plt.show()
plt.close()

png

We can group multiple galaxies at the same redshift into a Galaxies object, which is created from a list of individual galaxies.

In a strong lens, we ultimately group together a foreground lens galaxy and a background source galaxy. For now, we simply create a second galaxy and combine it with the original galaxy into a Galaxies object, to see how the light of multiple galaxies is represented.

extra_galaxy = al.Galaxy(
    redshift=0.5,
    bulge=al.lp.Sersic(
        centre=(0.2, 0.3),
        ell_comps=(0.0, 0.111111),
        intensity=1.0,
        effective_radius=1.0,
        sersic_index=2.5,
    ),
)

galaxies = al.Galaxies(galaxies=[galaxy, extra_galaxy])

The Galaxies object has similar methods to those for light profiles and individual galaxies.

For example, image_2d_from sums the images of all the galaxies.

image = galaxies.image_2d_from(grid=grid)

We can plot the combined image of all the galaxies, just like with other plotters.

aplt.plot_array(array=galaxies.image_2d_from(grid=grid), title="Image")

png

A subplot of each individual galaxy image can also be created.

aplt.subplot_galaxies(galaxies=galaxies, grid=grid)

png

Because galaxy light distributions often follow a log10 pattern, plotting in log10 space helps reveal details in the outskirts of the light profile.

This is especially helpful when visualizing how multiple galaxies overlap.

aplt.plot_array(array=galaxies.image_2d_from(grid=grid), title="Image", use_log10=True)

png

Units

Earlier, we mentioned that a galaxy's redshift allows us to convert between arcseconds and kiloparsecs.

A redshift measures how much a galaxy's light is stretched by the Universe's expansion. A higher redshift means the galaxy is further away, and its light has been stretched more. By knowing a galaxy's redshift, we can convert angular distances (like arcseconds) to physical distances (like kiloparsecs).

To perform this conversion, we use a cosmological model that describes the Universe's expansion. Below, we use the Planck15 cosmology, which is based on observations from the Planck satellite.

cosmology = al.cosmo.Planck15()

kpc_per_arcsec = cosmology.kpc_per_arcsec_from(redshift=galaxy.redshift)

print("Kiloparsecs per Arcsecond:")
print(kpc_per_arcsec)
Kiloparsecs per Arcsecond:
6.288247910157764

This kpc_per_arcsec can be used as a conversion factor between arcseconds and kiloparsecs when plotting images of galaxies.

We compute this value and plot the image, which by default is shown in units of arcseconds.

aplt.plot_array(array=galaxy.image_2d_from(grid=grid), title="Image")

png

Wrap Up

In this tutorial, you've learnt the basic quantities used to describe the galaxies that make up a strong lens, before we introduce any lensing calculations.

Let's summarise what we've covered:

  • Grids: A grid is a set of 2D $(y,x)$ coordinates that represent the positions where we measure the light of a galaxy.

  • Geometry: We showed how to shift, rotate, and convert grids to elliptical coordinates.

  • Light Profiles: Light profiles are analytic functions that describe how a galaxy's light is distributed in space. We used the Sersic profile to create images of galaxies.

  • Galaxies: Galaxies are collections of light profiles. We created galaxies with multiple light profiles, combined them into a Galaxies object, and visualized their images.

  • Units: By assuming redshifts for galaxies we can convert their quantities from arcseconds to physical units like kiloparsecs.

In the next tutorial, we'll introduce the mass of a galaxy and perform our first lensing calculation, whereby the light of a background source galaxy is deflected by the mass of a foreground lens galaxy.

Advanced Topics

The following advanced topics are not important for a new user learning the software for the first time. However, once you are an expert user, the following guides and concepts are important for doing accurate strong lens analysis, and thus may be things you want to commit to memory as future references.

Other Unit Conversion

Above, we used a redshift to convert between arcseconds and kiloparsecs. This is just one example of a unit conversion that can be performed using a galaxy's redshift.

There are many other unit conversions that can be performed, such as converting the units of a galaxy's image to what Astronomers call an AB magnitude system, which is a system used to measure the brightness of galaxies.

The autolens_workspace/*/guides/units module contains many examples of unit conversions and how to use them, but they will not be covered in the HowToLens tutorials.

Over Sampling

Over sampling is a numerical technique where the images of light profiles and galaxies are evaluated on a higher resolution grid than the image data to ensure the calculation is accurate.

For a new user, the details of over-sampling are not important, therefore just be aware that all calculations use an adaptive over sampling scheme with high accuracy across all use cases.

Once you are more experienced, you should read up on over-sampling in more detail via the autolens_workspace/*/guides/advanced/over_sampling.ipynb notebook.