Using Pyralysis

Using Pyralysis#

What this guide covers#

This guide provides a practical workflow for using Pyralysis for radio astronomy imaging (also called image synthesis / image formation). You’ll learn how to load data, simulate visibilities, grid, optimize (reconstruct) an image, and save results. Pyralysis often says reconstruction or optimization for the imaging step; see Glossary.

Tip

Best practice: Start with small datasets to get familiar with the workflow, then scale up to larger data using Dask. For pipeline workflows (SimulationPipeline, ImagerPipeline), optional local or cluster Dask setup is described in Pipelines and distributed Dask.

Note

Common mistake: Forgetting to set the correct cell size or image size for your data. Double-check these parameters!

—

Basic Workflow#

A typical workflow in Pyralysis involves:

  1. Loading Measurement Sets

  2. Simulating visibilities from an image (optional)

  3. Applying gridding techniques

  4. Imaging / reconstruction — run an optimization algorithm (or form a dirty image)

  5. Saving the final image

Note

Want the fastest way to generate a dataset? See Creating a dataset for the minimal example.

—

Loading Measurement Sets#

Pyralysis processes radio interferometric data stored in Measurement Set (MS) format. Under the hood, it leverages dask-ms for efficient handling of large datasets.

Example: Loading a Measurement Set

from pyralysis.io import DaskMS

# Define the path to your Measurement Set file
ms_file = "/path/to/data.ms"

# Load the dataset using Dask-MS (lazy loading enabled)
dataset = DaskMS(input_name=ms_file).read()

print("Loaded dataset structure:", dataset)

—

Simulating a Synthetic Dataset#

You can use Pyralysis to generate synthetic radio interferometric datasets for testing and benchmarking. This is useful for developing and validating imaging algorithms.

Example: Simulate an in-memory dataset (and optionally write a Measurement Set)

Note

The simulator returns a Dataset in memory. To choose the output .ms path, pass output_name to DaskMS and call write (not an argument on Simulator). See Creating a dataset and I/O Operations & Data Formats.

from pyralysis.io.antenna_config_io import AntennaConfigurationIo
from pyralysis.simulation import Simulator
from pyralysis.models.sky import PointSource
from pyralysis.io import DaskMS

# Load array configuration
interferometer = AntennaConfigurationIo(input_name="/path/to/array.cfg").read()
interferometer.configure_observation(
    min_frequency_hz=1e9, max_frequency_hz=1.1e9, frequency_step_hz=1e7,
    right_ascension="12:00:00", declination="45:00:00", integration_time=10, observation_time="1h"
)

# Define a simple sky model (point source)
source = PointSource(reference_intensity=1.0, sky_position="12:00:00 45:00:00")

# Run the simulation (Dataset stays in memory until you write)
sim = Simulator(interferometer=interferometer, sources=source)
dataset = sim.simulate(create_dataset=True)

# Write to disk — this is where you set the Measurement Set name/path
ms_io = DaskMS(output_name="/path/to/simulated.ms")
ms_io.write(dataset=dataset, tables="ALL")

Tip

See Simulation Framework for more advanced simulation options, including Gaussian sources, composite models, and noise injection.

—

Simulating Non-Parametric Sources#

Pyralysis allows simulating non-parametric sources, meaning we compute visibilities from an image without assuming any predefined model.

Methods for estimating visibilities from an image:

  • Nearest Neighbor Interpolation – Fastest but least accurate.

  • Bilinear Interpolation – Uses linear interpolation for smoother results.

  • Degridding (Convolutional Gridding with a Half-Kernel) – Most accurate method.

Example: Simulating Model Visibilities

from pyralysis.estimators import NearestNeighbor, BilinearInterpolation, Degridding
from pyralysis.convolution import PSWF1

# Use Degridding (Convolutional Gridding)
kernel = PSWF1(size=3, cellsize=0.003, oversampling_factor=3)
mv_degridding = Degridding(input_data=dataset, image=image, cellsize=0.003, ckernel_object=kernel)

# Compute visibilities
mv_degridding.transform()

—

Gridding Visibilities#

Gridding is essential for transforming irregularly sampled visibilities into a uniform Fourier grid. Pyralysis provides DirtyMapper to compute a dirty image and dirty beam.

Example: Applying Gridding

from pyralysis.transformers import DirtyMapper
from pyralysis.convolution import PSWF1

# Define a convolution kernel
kernel = PSWF1(size=3, cellsize=0.003, oversampling_factor=3)

# Create a DirtyMapper object for gridding
dirty_mapper = DirtyMapper(input_data=dataset, imsize=512, cellsize=0.003, padding_factor=1.2, ckernel_object=kernel)

# Compute dirty image and dirty beam
dirty_images, dirty_beam = dirty_mapper.transform()

Both dirty_images and dirty_beam are returned as Image objects.

🔹 For more details on gridding, see Gridding Techniques in Pyralysis. For image processing and I/O operations, see Image Processing and I/O Operations & Data Formats.

—

Performing Image Reconstruction#

To reconstruct an image, Pyralysis uses gradient-based optimization techniques.

Example: Running an Optimizer

from pyralysis.optimization.optimizer import HagerZhang
from pyralysis.optimization.terms import ChiSquared, L1Norm
from pyralysis.optimization import ObjectiveFunction

# Define the objective function
term_list = [ChiSquared(model_visibility=mv, normalize=True), L1Norm(penalization_factor=0.005)]
of = ObjectiveFunction(term_list=term_list, image=image, persist_gradient=True)

# Run the optimizer
optim = HagerZhang(image=image, objective_function=of, max_iter=100)
reconstructed_image = optim.optimize()

🔹 For details on optimizers, line searchers, and seeders, see Optimizers, Line Searchers, and Step Size Seeders (overview: Optimization in Pyralysis).

—

Phase Center Shifting#

Pyralysis provides the Shifter transformer for changing the effective phase center of observations. This is essential for:

  • Multi-field imaging where different fields have different phase centers

  • Phase center corrections for proper source positioning

  • Wide-field imaging with non-coplanar baselines

  • Source positioning in visibility space

Example: Shifting Phase Center

from pyralysis.transformers import Shifter
from astropy.coordinates import SkyCoord
import astropy.units as u

# Shift to a new phase center
new_center = SkyCoord(ra=12.5*u.deg, dec=-45.2*u.deg)
shifter = Shifter(phase_center=new_center)
shifter.transform()

# Shift relative to a specific reference
ref_center = SkyCoord(ra=12.0*u.deg, dec=-45.0*u.deg)
shifter = Shifter(phase_center=new_center, reference_phase_center=ref_center)
shifter.transform()

The phase shift is applied using the full UVW coordinate system:

\[V'(u,v,w) = V(u,v,w) \cdot e^{-2\pi i (u l_0 + v m_0 + w (n_0 - 1))}\]

where \((l_0, m_0, n_0)\) are the direction cosines of the phase center offset.

Tip

The Shifter automatically handles visibility-like columns (DATA, MODEL_DATA, CORRECTED_DATA, RESIDUAL_DATA, and other *_DATA columns) and updates both the visibility data and field phase center information.

—

Saving Reconstructed Images#

Pyralysis provides I/O handlers to save reconstructed images.

Example: Saving Images in FITS or Zarr Format

from pyralysis.io import FITS, ZarrArray

ioh_fits = FITS()  # FITS format
ioh_zarr = ZarrArray()  # Zarr format

# Use the chosen I/O handler when running the optimizer
optim = HagerZhang(image=image, objective_function=of, max_iter=100, io_handler=ioh_fits)

🔹 For details on I/O handlers, see pyralysis.

—

Next Steps#

Now that you’ve seen the basic usage of Pyralysis, explore the advanced topics:

See also#

Quickstart Guide | Optimization in Pyralysis