karabo.simulation.signal package

Submodules

karabo.simulation.signal.base_segmentation module

Base segmentation class.

class BaseSegmentation

Bases: ABC

Base segmentation class.

abstract segment(image: Image3D) SegmentationOutput

Segments on a 2D or 3D image.

karabo.simulation.signal.base_signal module

Base signal class.

class BaseSignal

Bases: Generic[T], ABC

Base signal class.

abstract simulate() list[T]

Simulate a signal to get a 2D/3D image output.

karabo.simulation.signal.eor_profile module

EoR profile simulation.

class EoRProfile

Bases: object

EoR profile simulation.

Examples

>>> rnd = np.random.rand(200)
>>> sum_rnd = np.cumsum(rnd)
>>> sum_rnd /= sum_rnd[-1]
>>> eor = EoRProfile.simulate(x_hi=sum_rnd)
>>> EoRProfile.plot(eor)
>>> plt.show()
delta_m = 0

Density contrast of matter.

frequency_21cm = 1420405751.768

Frequency of the 21cm signal [Hz].

hubble_constant = 0.7

Hubble constant.

hubble_param = 70.0

Hubble parameter [km/s/Mpc].

omega_b = 0.05

Baryon density parameter.

omega_m = 0.31

Matter density parameter.

classmethod plot(x_hi: Annotated[ndarray[Any, dtype[float64]], Literal['N']] | None = None, profile: Annotated[ndarray[Any, dtype[float64]], Literal['N', 2]] | None = None) Figure

Plot the fluctuation profile of the 21cm signal.

Either the x_hi parameter needs to be given or the calculated profile.

Parameters

x_hiOptional[Annotated[npt.NDArray[np.float_], Literal[“N”]]], optional

Neutral hydrogen fraction, by default None.

profileOptional[EoRProfileT], optional

An optional profile to be plotted. If not given, a default EoR profile will be plotted. By default None.

Returns

Figure

The plotted figure of the 21cm signal.

classmethod simulate(x_hi: Annotated[ndarray[Any, dtype[float64]], Literal['N']], dv_r_over_dr: float = 0, f_range: tuple[float, float] = (1000000.0, 200000000.0)) Annotated[ndarray[Any, dtype[float64]], Literal['N', 2]]

Calculate the approximate evolution of fluctuations in the 21cm brightness.

The number of points is determined by the resolution of the x_hi array.

Implemented per https://arxiv.org/pdf/1602.02351.pdf, equation (1)

Parameters

x_hiAnnotated[npt.NDArray[np.float_], Literal[“N”]]

Neutral hydrogen fraction.

dv_r_over_drfloat, optional

???. By default 0

f_rangetuple[float, float], optional

Frequency range to plot in [Hz]. by default (2, 200e6)

Returns

EoRProfileT

An array of the shape (floor((f_end - f_start) / step_size), 2), containing the frequency in the first column and the corresponding EoR profile in the second.

t_gamma = 2.73

Cosmic microwave background temperature [Kelvin].

t_s = 0.068

Spin temperature [Kelvin].

karabo.simulation.signal.galactic_foreground module

Galactic Foreground signal catalogue wrapper.

class SignalGalacticForeground(centre: SkyCoord, fov: Angle, redshifts: list[float], grid_size: tuple[Annotated[int, Literal['X']], Annotated[int, Literal['Y']]], gleam_file_path: Path | None = None)

Bases: BaseSignal[Image2D]

Galactic Foreground signal catalogue wrapper.

Examples

>>> from karabo.simulation.signal.plotting import SignalPlotting
>>> available_redshifts = SignalGalacticForeground.available_redshifts()
>>> cent = SkyCoord(ra=10 * units.degree, dec=20 * units.degree, frame="icrs")
>>> gf = SignalGalacticForeground(
...    cent,
...    redshifts=available_redshifts[:3],
...    grid_size=(30, 30),
...    fov=Angle([20, 20], unit=units.degree),
... )
>>> imgs = gf.simulate()
>>> SignalPlotting.general_img(imgs[0], "Galactic foreground", log_bar=True)
>>> # Map plotting
>>> col = SignalGalacticForeground._flux_column(7.6)
>>> data = (
...     gf.gleam_catalogue[gf.gleam_catalogue[col] >= 0]
...         .to_pandas()
...         .sort_values(SignalGalacticForeground.RA_COLUMN)
... )
>>> SignalPlotting.general_polar_plot(
...     data[SignalGalacticForeground.RA_COLUMN],
...     data[SignalGalacticForeground.DEC_COLUMN],
...     data[col],
...     title="Galactic Foreground",
...     log_bar=True,
... )
DEC_COLUMN: Final[str] = 'DEJ2000'
RA_COLUMN: Final[str] = 'RAJ2000'
classmethod available_redshifts(gleam_file_path: Path | None = None, gleam_catalogue: Table | None = None) list[float]

Get a list from available redshift values (as float).

If the gleam catalogue is given, it will take precedents over the path. Else, if no path is given, the default catalogue will bi used.

Parameters

gleam_file_pathOptional[Path], optional

Path to the gleam catalogue, by default None

gleam_catalogueOptional[Table], optional

Gleam catalogue data, by default None

Returns

list[float]

List of available redshift values.

simulate() list[Image2D]

Simulate a signal to get a 2D image output.

karabo.simulation.signal.helpers module

Helpers for generating different types of signals.

filter_dataframe_radec(df: DataFrame, centre: SkyCoord, fov: Angle, ra_column: str, dec_column: str, wrap_offset: bool = False) tuple[DataFrame, Angle, Angle]

Filter a dataframe with RA-DEC coordinates.

Expects the Center to be at (0°, 0°) with the RA axis having a range of (-180°, 180°) and the DEC axis (-90°, 90°).

Parameters

dfpd.DataFrame

The dataframe that is to be filtered.

centreSkyCoord

The center position of the filter.

fovAngle

Field of View in the (RA, DEC) axis.

ra_columnstr

Name of the RA column in the dataframe.

dec_columnstr

Name of the DEC column in the dataframe.

wrap_offsetbool

If we reach a wrap point, if all the coordinates should be wrapped back into the positive space.

Returns

tuple[pd.DataFrame, Angle, Angle]

The filtered dataframe (copy), lower left corner, upper right corner.

interpolate_image(image: Annotated[ndarray[Any, dtype[float64]], Literal['X', 'Y']], new_size: tuple[int, int], mode: Literal['reflect', 'grid-mirror', 'constant', 'grid-constant', 'nearest', 'mirror', 'grid-wrap', 'wrap'] = 'constant') Annotated[ndarray[Any, dtype[float64]], Literal['N', 'M']]

Interpolate an image from it’s original size to new_size.

Parameters

imageAnnotated[npt.NDArray[np.float_], Literal[“X”, “Y”]]

The image that is to be interpolated.

new_sizetuple[int, int]

The desired output size of the image

modeSciPyInterpolationModes, optional

The interpolation mode, by default “constant”.

Returns

Annotated[npt.NDArray[np.float_], Literal[“N”, “M”]]

The interpolated input image with it’s new dimensions.

map_radec_datapoints_to_grid(data: DataFrame, grid_size: tuple[Annotated[int, Literal['X']], Annotated[int, Literal['Y']]], ra_column: str, dec_column: str, intensity_column: str, par_count: int = 5) Annotated[ndarray[Any, dtype[float64]], Literal['X', 'Y']]

Map the given data points with a destination to source mapping.

For each pixel in the destination grid, the equivalent degree range in the source will be summed together and set in the destination grid.

Parameters

datapd.DataFrame,

The data that is to be plotted onto the grid.

grid_sizetuple[Annotated[int, Literal[“X”]], Annotated[int, Literal[“Y”]]]

Size of the output grid.

ra_columnstr

Name of the RA column in the dataframe.

dec_columnstr

Name of the DEC column in the dataframe.

intensity_columnstr

Name of the column to use for the intensities.

par_countint, optional

How many processes to use to map the data points. By default 5.

Returns

Annotated[npt.NDArray[np.float_], Literal[“X”, “Y”]]

A 2D numpy array representing an image with the dimensions of the grid_size parameter.

karabo.simulation.signal.plotting module

Signal plotting helpers.

class SegmentationPlotting

Bases: object

Plotting utilities for the segmentation.

classmethod seg_u_net_plotting(segmented: SegmentationOutput, bin_count: int = 20) Figure

Plot the first slice of the segU-net cube.

Parameters

segmentedSegmentationOutput

Output of the segmentation

bin_countint

The number output bins in the histogram plot

Returns

Figure

Plot of the SegU-Net segmentation

Raises

KaraboError

If the input ‘xhi_seg_err’ is None

classmethod superpixel_plotting(segmented: SegmentationOutput, signal_image: Image3D, log_sky: bool = False) Figure

Plot the first slice of the superpixel cube.

Parameters

segmentedSegmentationOutput

output of the segmentation

signal_imageImage3D

Image cube

log_skybool, optional

If the color bar of the sky plot should have a symmetric log norm applied.

Returns

Figure

Plot of the Superpixel segmentation

Raises

KaraboError

If the input ‘xhii_stitch’ is None

class SignalPlotting

Bases: object

Signal plotting helpers.

classmethod brightness_temperature(data: BaseImage, z_layer: int = 0) Figure

Plot the brightness temperature of a 2D image.

Parameters

dataBaseImage

The image to be plotted.

z_layerint, optional

The Z layer to be used, when a Image3D is used.

Returns

Figure

Figure of the plotted image

classmethod general_img(img: Image2D, title: str, tick_count: int = 5, x_label: str = 'RA [°]', y_label: str = 'DEC [°]', bar_label: str = 'Temperature [K]', log_bar: bool = False) Figure

Plot a general image with a temperature.

Parameters

imgImage2D

The image to be plotted.

titlestr

Title to be shown in the figure.

tick_countint, optional

The count of ticks to show along each axis, by default 5

x_labelstr, optional

Label to be plotted along the X-axis.

y_labelstr, optional

Label to be plotted along the Y-axis.

bar_label: str, optional

Label for the color bar.

log_barbool, optional

If the color bar should have a symmetric log norm applied.

Returns

Figure

The resulting plot figure.

classmethod general_polar_plot(ra_series: Annotated[ndarray[Any, dtype[float64]], Literal['N']], dec_series: Annotated[ndarray[Any, dtype[float64]], Literal['N']], intensities_series: Annotated[ndarray[Any, dtype[float64]], Literal['N']], title: str, bar_label: str = 'Temperature [K]', log_bar: bool = False) Figure

Plot a RA/DEC data in a polar plot with the intensity representing the color.

Parameters

ra_seriesAnnotated[npt.NDArray[np.float_], Literal[“N”]]

RA coordinates in degrees.

dec_seriesAnnotated[npt.NDArray[np.float_], Literal[“N”]]

DEC coordinates in degrees.

intensities_seriesAnnotated[npt.NDArray[np.float_], Literal[“N”]]

Intensities in Kelvin.

titlestr

Title to be shown in the figure.

bar_label: str, optional

Label for the color bar.

log_barbool, optional

If the color bar should have a symmetric log norm applied.

Returns

Figure

The resulting plot figure.

classmethod power_spectrum(data: BaseImage | SegmentationOutput | list[BaseImage | SegmentationOutput], kbins: int = 15, z_layer: int = 0) Figure

Plot the power spectrum the 21cm signal.

Parameters

dataUnion[BaseImage, SegmentationOutput]

Either a BaseImage or a SegmentationOutput

kbinsint, optional

Count of bins for the spectrum plot, by default 15

z_layerint

If an Image3D is passed, then this layer will be used to plot the image. By default 0

Returns

Figure

The generated plot figure.

classmethod power_spectrum_image_vs_xfrac_dens(image: BaseImage | SegmentationOutput, xfrac_dens: XFracDensFilePair, kbins: int = 15, z_layer: int = 0) Figure

Plot the power spectrum by using an image and the original xfrac/dens file pair.

Parameters

imageUnion[BaseImage, SegmentationOutput]

Either a BaseImage or a SegmentationOutput

xfrac_densXFracDensFilePair

The xfrac and dens file pair.

kbinsint, optional

Count of bins for the spectrum plot, by default 15

z_layerint

If an Image3D is passed, then this layer will be used to plot the image. By default 0

Returns

Figure

The generated plot figure.

classmethod power_spectrum_xfrac_dens(data: XFracDensFilePair | list[XFracDensFilePair], kbins: int = 15) Figure

Plot the power spectrum the 21cm signal using xfrac and dens-files.

Parameters

dataUnion[XFracDensFilePair, list[XFracDensFilePair]]

The xfrac and dens file pair or a list thereof.

kbinsint, optional

Count of bins for the spectrum plot, by default 15

Returns

Figure

The generated plot figure.

classmethod xfrac_dens(data: XFracDensFilePair) Figure

Plot the xfrac and dens files.

Parameters

dataXFracDensFilePair

The xfrac and dens file pair.

Returns

Figure

The figure that was plotted.

karabo.simulation.signal.seg_u_net_segmentation module

Segmentation with SegU-net.

class FixedSegUNet(tta: Literal[0, 1, 2] = 1, verbose: bool = False)

Bases: segunet21cm

Fixes the bug such that the Segmentation can run.

Can be removed once Tools21cm is fixed.

class SegUNetSegmentation(max_baseline: float = 70.0, tta: Literal[0, 1, 2] = 2)

Bases: BaseSegmentation

SegU-net based segmentation. Using its own init for t2c (Tools21cm).

Examples

>>> from karabo.simulation.signal.plotting import SegmentationPlotting
>>> from karabo.simulation.signal.signal_21_cm import Signal21cm
>>> z = Signal21cm.get_xfrac_dens_file(z=7.059, box_dims=244 / 0.7)
>>> sig = Signal21cm([z])
>>> signal_images = sig.simulate()
>>> seg = SegUNetSegmentation(max_baseline=70.0, tta=2)
>>> segmented = seg.segment(signal_images[0])
>>> SegmentationPlotting.seg_u_net_plotting(segmented=segmented)
segment(image: Image3D) SegmentationOutput

SegU-net based segmentation. Using its own init for t2c (Tools21cm).

Parameters

imageImage3D

The constructed simulation

Returns

SegmentationOutput

SegU-net cube

karabo.simulation.signal.signal_21_cm module

21cm Signal simulation.

class Signal21cm(files: list[XFracDensFilePair])

Bases: BaseSignal[Image3D]

21cm Signal simulation wrapper.

Examples

>>> from karabo.simulation.signal.plotting import SignalPlotting
>>> redshifts = Signal21cm.available_redshifts()
>>> z1 = Signal21cm.get_xfrac_dens_file(z=redshifts[0], box_dims=244 / 0.7)
>>> z2 = Signal21cm.get_xfrac_dens_file(z=redshifts[1], box_dims=244 / 0.7)
>>> sig = Signal21cm([z1, z2])
>>> signal_images = sig.simulate()
>>> fig = SignalPlotting.brightness_temperature(signal_images[0])
>>> fig.savefig("brightness_temperature.png")
>>> from karabo.simulation.signal.plotting import SignalPlotting
>>> z = Signal21cm.randomized_lightcones(200, 8)
>>> SignalPlotting.brightness_temperature(z, z_layer=100)
astro_gmell_base_url = 'https://ttt.astro.su.se/~gmell/244Mpc'
classmethod available_redshifts() list[float]

Get all available redshifts.

Returns

list[float]

List of all available redshifts.

classmethod available_redshifts_dens() list[float]

Get all available redshifts for dens files.

Returns

list[float]

List of all available redshifts for dens files.

classmethod available_redshifts_xfrac() list[float]

Get all available redshifts for xfrac files.

Returns

list[float]

List of all available redshifts for xfrac files.

static default_r_hii(redshift: float) float

Light cone HII region size calculation function (default implementation).

Parameters

redshiftfloat

Redshift, to determine the radius for.

Returns

float

Light cone radius.

dens_url = 'densities/nc250/coarser_densities'
static get_xfrac_dens_file(z: float, box_dims: float) XFracDensFilePair

Get the xfrac and dens files from the server.

They are downloaded and cached on the first access.

Parameters

zfloat

Redshift value.

box_dimsfloat

Box dimensions used for these files.

Returns

XFracDensFilePair

A tuple of xfrac and dens files.

classmethod randomized_lightcones(n_cells: int, z: float, r_hii: Callable[[float], float] | None = None) Image3D

Generate an image with randomized light cones.

Parameters

n_cellsint

The count of cells to produce.

zfloat

The redshift value for this image.

r_hiiCallable[[float], float], optional

Radius function of the HII region. By default None, resulting in the execution of Signal21cm.default_r_hii

Notes

Implementation according to https://tools21cm.readthedocs.io/examples/lightcone.html

Returns

Image3D

The generated cube with multiple light cones.

simulate() list[Image3D]

Simulate the 21cm signal as a 3D intensity cube.

Returns

list[Image3D]

A list of 3D image cubes, based on the self.files list of provided xfrac and dens files. The pixel values are in Kelvin.

Raises

KaraboError

If a pair of xfrac and dens files do not have the same redshift values.

xfrac_url = '244Mpc_f2_0_250'

karabo.simulation.signal.superimpose module

Superimpose two or more signals.

class Superimpose

Bases: object

Superimpose two or more signals.

classmethod combine(*signals: BaseImage) BaseImage

Superimpose two or more signals int a single signal.

Superimposing is done by adding each signal to the previous one. To combine a Image3D with a Image2D, the Image2D will be added to every layer of the Image3D.

If only one signal is passed, it gets returned without any further processing.

Parameters

signalsBaseImage

The signals that are to be combined.

Returns

BaseImage

If the input includes at least one Image3D, the return type is a Image3D, else Image2D.

Raises

KaraboError

When an empty signal list is passed in.

karabo.simulation.signal.superpixel_segmentation module

Segmentation with Superpixel.

class SuperpixelSegmentation(max_baseline: float = 70.0, n_segments: int = 1000, max_iter: int = 5)

Bases: BaseSegmentation

Superpixel based segmentation.

Examples

>>> from karabo.simulation.signal.plotting import SegmentationPlotting
>>> from karabo.simulation.signal.signal_21_cm import Signal21cm
>>> z1 = Signal21cm.get_xfrac_dens_file(z=7.059, box_dims=244 / 0.7)
>>> sig = Signal21cm([z1])
>>> seg = SuperpixelSegmentation(max_baseline=70.0,  max_iter=5, n_segments=1000)
>>> signal_images = sig.simulate()
>>> segmented = seg.segment(signal_images[0])
>>> SegmentationPlotting.superpixel_plotting(segmented, signal_images[0])
segment(image: Image3D) SegmentationOutput

Superpixel based segmentation.

Parameters

imageImage3D

The constructed simulation

Returns

SegmentationOutput

Superpixel cube

karabo.simulation.signal.synchroton_signal module

Synchroton signal catalogue wrapper.

class SignalSynchroton(centre: SkyCoord, fov: Angle, grid_size: tuple[Annotated[int, Literal['X']], Annotated[int, Literal['Y']]], diffuse_emission_path: Path | None = None, data_table_id: int = 1, mask_table_id: int = 2)

Bases: BaseSignal[Image2D]

Synchroton signal.

Examples

>>> from karabo.simulation.signal.plotting import SignalPlotting
>>> from astropy import units
>>> cent = SkyCoord(ra=10 * units.degree, dec=10 * units.degree, frame="icrs")
>>> sync_sig = SignalSynchroton(
...     centre=cent,
...     fov=Angle([20, 20], unit=units.degree),
...     grid_size=(100, 100),
... )
>>> imgs = sync_sig.simulate()
>>> SignalPlotting.general_img(imgs[0], "Synchroton")
>>> # Map overview plot
>>> SignalPlotting.general_polar_plot(
...     sync_sig.data["RA"],
...     sync_sig.data["DEC"],
...     sync_sig.data["intensity"],
...     title="Synchroton map",
...     log_bar=True,
... )
MAX_GRID_X = 190
MAX_GRID_Y = 190
simulate() list[Image2D]

Fetch the synchroton signal.

Returns

list[Image2D]

The synchroton image.

karabo.simulation.signal.typing module

General typings for the signal package.

class BaseImage(data: Annotated[ndarray[Any, dtype[float64]], Literal['X', 'Y']], x_label: Annotated[ndarray[Any, dtype[float64]], Literal['X']], y_label: Annotated[ndarray[Any, dtype[float64]], Literal['Y']], redshift: float, box_dims: float)

Bases: object

A general image, meant to be subclassed.

box_dims: float

Box dimensions used to create the image. Negative number if not provided.

data: Annotated[ndarray[Any, dtype[float64]], Literal['X', 'Y']]

Image data in a 2D array.

redshift: float

Redshift value. Negative number if not provided.

x_label: Annotated[ndarray[Any, dtype[float64]], Literal['X']]

X-labels.

y_label: Annotated[ndarray[Any, dtype[float64]], Literal['Y']]

Y-labels.

class Image2D(data: Annotated[ndarray[Any, dtype[float64]], Literal['X', 'Y']], x_label: Annotated[ndarray[Any, dtype[float64]], Literal['X']], y_label: Annotated[ndarray[Any, dtype[float64]], Literal['Y']], redshift: float, box_dims: float)

Bases: BaseImage

A 2D image.

class Image2DOriented(data: Annotated[ndarray[Any, dtype[float64]], Literal['X', 'Y']], x_label: Annotated[ndarray[Any, dtype[float64]], Literal['X']], y_label: Annotated[ndarray[Any, dtype[float64]], Literal['Y']], redshift: float, box_dims: float, sky_model: SkyModel)

Bases: Image2D

A 2D image combined with a sky model.

sky_model: SkyModel

Sky Model for the orientation.

class Image3D(data: Annotated[ndarray[Any, dtype[float64]], Literal['Z', 'X', 'Y']], x_label: Annotated[ndarray[Any, dtype[float64]], Literal['X']], y_label: Annotated[ndarray[Any, dtype[float64]], Literal['Y']], redshift: float, box_dims: float, z_label: Annotated[ndarray[Any, dtype[float64]], Literal['Z']])

Bases: BaseImage

A 3D cube of images along the z-axis.

data: Annotated[ndarray[Any, dtype[float64]], Literal['Z', 'X', 'Y']]

Image data in a 3D cube.

z_label: Annotated[ndarray[Any, dtype[float64]], Literal['Z']]

Z-labels.

class SegmentationOutput(image: Image3D, xhii_stitch: npt.NDArray[np.bool_] | None, mask_xhi: npt.NDArray[np.bool_], dt_smooth: npt.NDArray[np.float_], xhi_seg_err: npt.NDArray[np.float_] | None)

Bases: NamedTuple

Output of the segmentation.

dt_smooth: ndarray[Any, dtype[float64]]

Alias for field number 3

image: Image3D

Alias for field number 0

mask_xhi: ndarray[Any, dtype[bool_]]

Alias for field number 2

xhi_seg_err: ndarray[Any, dtype[float64]] | None

Alias for field number 4

xhii_stitch: ndarray[Any, dtype[bool_]] | None

Alias for field number 1

class XFracDensFilePair(xfrac_path: Path, dens_path: Path, box_dims: float)

Bases: NamedTuple

A pair of matching XFrac and dens files (Same Redshift value).

box_dims: float

Length of the volume along each direction in [Mpc].

dens_path: Path

Path to the dens file.

load() XFracDensLoaded

Load the files into memory.

Returns

XFracDensLoaded

The loaded files.

xfrac_path: Path

Path to the x-frac file.

class XFracDensLoaded(x_file: t2c.XfracFile, d_file: t2c.DensityFile, x_frac: FloatArrayNxNxN, dens: FloatArrayNxNxN, z: float, box_dims: float)

Bases: NamedTuple

The Xfrac and dens files loaded into memory.

box_dims: float

Length of the volume along each direction in [Mpc].

d_file: DensityFile

Loaded density file object.

dens: Annotated[ndarray[Any, dtype[float64]], Literal['N']]

Density values in a 3D cub N*N*N.

x_file: XfracFile

Loaded xfrac file object.

x_frac: Annotated[ndarray[Any, dtype[float64]], Literal['N']]

xfrac values in a 3D cube N*N*N.

xy_dims() tuple[Annotated[ndarray[Any, dtype[float64]], Literal['N']], Annotated[ndarray[Any, dtype[float64]], Literal['N']]]

Get the x, y dimensions of the loaded.

Returns

tuple[FloatArrayN, FloatArrayN]

A tuple containing the x and y labels.

z: float

Redshift value

Module contents