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.
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.
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
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.
- 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
- 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