Full Class List

BayesEoRParser

class bayeseor.params.BayesEoRParser(*args, parse_as_dict=False, **kwargs)

Command line argument parser class for all BayesEoR analysis parameters.

For more information on the command line syntax, please see

python run-analysis.py --help
General Parameters
config

Path to a yaml configuration file containing any command line arguments below parsable by jsonargparse.ArgumentParser. Please note command line arguments use dashes ‘-’ but the configuration yaml file requires underscores ‘_’. For example, to specify the array directory prefix on the command line, the command line argument is –array-dir-prefix. To specify the array directory prefix in the configuration yaml file, use array_dir_prefix instead. For an example, please see the provided configuration yaml file BayesEoR/example-config.yaml.

Type:

str, optional

quiet

Quiet mode where only a few statements are printed. Defaults to False.

Type:

bool, optional

clobber

Overwrite files if they exist. This includes the matrix stack, any existing preprocessed data/noise vectors, and instrument model files. Defaults to False.

Type:

bool, optional

Compute Parameters
use_gpu

Use GPUs (True) or CPUs (False). Defaults to True.

Type:

bool, optional

run

Run a full power spectrum analysis. To build the matrix stack only, run can be omitted or set to False in the configuration yaml.

Type:

bool, optional

Matrix Parameters
array_dir_prefix

Directory for matrix storage. Defaults to ‘./array-storage’.

Type:

str, optional

use_sparse_matrices

Use sparse matrices (True) to reduce storage requirements or dense matrices (False). Defaults to True.

Type:

bool, optional

build_Finv_and_Fprime

If True, construct Finv and Fprime independently and write both matrices to disk when building the matrix stack. Otherwise (default), construct the matrix product Finv_Fprime in place from the dense matrices comprising Finv and Fprime to minimize the memory and time required to build the matrix stack. In this case, only the matrix product Finv_Fprime is written to disk.

Type:

bool, optional

Prior Parameters
priors

Power spectrum prior range [min, max] for each k bin.

Type:

list of list

log_priors

Assume priors on power spectrum coefficients are in log_10 units (True) or linear units (False). Defaults to True.

Type:

bool, optional

uprior_bins

Array indices of k-bins using a uniform prior. All other bins use the default log-uniform prior. Follows python slicing syntax. Can pass a range via ‘1:4’ (non-inclusive high end), a list of indices via ‘1,4,6’ (no spaces between commas), a single index ‘3’ or ‘-3’, or ‘all’. Defaults to an empty string (all k-bins use log-uniform priors).

Type:

str, optional

inverse_LW_power

Prior on the inverse power of the large spectral scale model (LSSM) coefficients. A large value, 1e16, constrains the LSSM coefficients to be zero. A small value, 1e-16 (default), leaves the LSSM coefficients unconstrained.

Type:

float, optional

use_LWM_Gaussian_prior

Use a Gaussian prior on the large spectral scale model (True, NOT CURRENTLY IMPLEMENTED). Otherwise, use a uniform prior (False, default)

Type:

bool, optional

Sampler Parameters
output_dir

Directory for sampler output. Defaults to ‘./chains/’.

Type:

str, optional

file_root

If None (default), start a new analysis. Otherwise, resume analysis from file_root.

Type:

str, optional

use_Multinest

Use MultiNest (True) or Polychord (False) as the sampler. Using Polychord is advised for large parameter spaces. Defaults to True (use MultiNest).

Type:

bool, optional

Frequency Parameters
nf

Number of frequency channels. Required if data_path points to a preprocessed data vector with a ‘.npy’ suffix. Otherwise, nf sets the number of frequencies to keep starting from freq_idx_min or freq_min, or around freq_center. Defaults to None (keep all frequencies).

Type:

int, optional

df

Frequency channel width in hertz. Required if data_path points to a preprocessed numpy-compatible file. Otherwise, if None (default), defaults to the frequency channel width in the input pyuvdata-compatible visibilities. Defaults to None.

Type:

float, optional

freq_idx_min

Minimum frequency channel index to keep in the data vector. Used only if data_path points to a pyuvdata-compatible visibility file. Defaults to None (keep all frequencies).

Type:

int, optional

freq_min

Minimum frequency in hertz. If data_path points to a pyuvdata-compatible visibility file, freq_min sets the minimum frequency kept in the data vector. All frequencies greater than or equal to freq_min will be kept, unless nf is specified. If None (default), all frequencies are kept. If data_path points to a preprocessed data vector with a ‘.npy’ suffix, one of freq_min or freq_center is required.

Type:

float, optional

freq_center

Central frequency in Hertz. If data_path points to a pyuvdata-compatible visibility file, nf is also required to determine the number of frequencies kept around freq_center in the data vector. If None (default), all frequencies are kept. If data_path points to a preprocessed data vector with a ‘.npy’ suffix, one of freq_min or freq_center is required.

Type:

float, optional

neta

Number of line-of-sight Fourier modes. Defaults to nf.

Type:

int, optional

nq

Number of large spectral scale model quadratic basis vectors. If beta is not None, the quadratic basis vectors are replaced by power law basis vectors according to the spectral indices in beta. Defaults to 0.

Type:

int, optional

beta

Brightness temperature power law spectral index/indices used in the large spectral scale model. Can be a single spectral index, ‘[2.63]’, or multiple spectral indices can be passed, ‘[2.63,2.82]’, to use multiple power law spectral basis vectors. Do not put spaces after commas if using multiple spectral indices. Defaults to [2.63, 2.82].

Type:

list of float, optional

fit_for_spectral_model_parameters

Fit for the optimal large spectral scale model spectral indices. Defaults to False.

Type:

bool, optional

pl_min

Minimum brightness temperature spectral index when fitting for the optimal large spectral scale model spectral indices. Defaults to None.

Type:

float, optional

pl_max

Maximum brightness temperature spectral index when fitting for the optimal large spectral scale model spectral indices. Defaults to None.

Type:

float, optional

pl_grid_spacing

Grid spacing for the power law spectral index axis when fitting for the optimal large spectral scale model spectral indices. Defaults to None.

Type:

float, optional

Time Parameters
nt

Number of times. Required if data_path points to a preprocessed data vector with a ‘.npy’ suffix. Otherwise, sets the number of times to keep starting from jd_idx_min or jd_min, or around jd_center. Defaults to None (keep all times).

Type:

int, optional

dt

Integration time in seconds. Required if data_path points to a preprocessed data vector with a ‘.npy’ suffix. Otherwise, if None (default), defaults to the integration time in the input pyuvdata-compatible visibilities. Defaults to None.

Type:

float, optional

jd_idx_min

Minimum time index to keep in the data vector. Used only if data_path points to a pyuvdata-compatible visibility file. Defaults to None (keep all times). Defaults to None (keep all times).

Type:

int, optional

jd_min

Minimum time as a Julian date. If data_path points to a pyuvdata-compatible visibility file, jd_min sets the minimum time kept in the data vector. All times greater than or equal to jd_min will be kept, unless nt is specified. If None (default), all times are kept. If data_path points to a preprocessed data vector with a ‘.npy’ suffix, one of jd_min or jd_center is required.

Type:

float, optional

jd_center

Central time as a Julian date. If data_path points to a pyuvdata-compatible visibility file, nt is also required to determine the number of times kept around jd_center in the data vector. If None (default), all times are kept. If data_path points to a preprocessed data vector with a ‘.npy’ suffix, one of jd_min or jd_center is required.

Type:

float, optional

Model Image Parameters
nside

HEALPix resolution parameter. Sets the resolution of the sky model. Note, the HEALPix resolution must be chosen such that there are two HEALPix pixels per minimum fringe wavelength from the model uv-plane to satisfy the Nyquist-Shannon sampling theorem.

Type:

int

fov_ra_eor

Field of view of the right ascension axis of the EoR sky model in degrees.

Type:

float

fov_dec_eor

Field of view of the declination axis of the EoR sky model in degrees. Defaults to fov_ra_eor.

Type:

float, optional

fov_ra_fg

Field of view of the right ascension axis of the foreground sky model in degrees. Defaults to fov_ra_eor.

Type:

float, optional

fov_dec_fg

Field of view of the declination axis of the foreground sky model in degrees. Defaults to fov_ra_fg or fov_dec_eor if fov_ra_fg is not defined.

Type:

float, optional

simple_za_filter

Filter pixels in the sky model by zenith angle only (True, default). Otherwise, filter pixels in a rectangular region set by the field of view values along the RA and Dec axes (False). It is suggested to set simple_za_filter to True (see issue #11).

Type:

bool, optional

Model uv-Plane Parameters
nu

Number of pixels on the u-axis of the model uv-plane for the EoR model.

Type:

int

nv

Number of pixels on the v-axis of the model uv-plane for the EoR model. Defaults to nu.

Type:

int, optional

nu_fg

Number of pixels on the u-axis of the model uv-plane for the foreground model. Defaults to nu.

Type:

int, optional

nv_fg

Number of pixels on the v-axis of the model uv-plane for the foreground model. Defaults to nu_fg or nv if nu_fg is not defined.

Type:

int, optional

fit_for_monopole

Include (True) or exclude (False) the (u, v) == (0, 0) pixel in the model uv-plane. Defaults to False.

Type:

bool, optional

Noise Model Parameters
sigma

Standard deviation of the visibility noise in mK sr. Required if calc_noise is False and data_path points to a pyuvdata-compatible visibility file or noise_data_path is None and data_path points to a preprocessed numpy-compatible visibility vector. Defaults to None.

Type:

float, optional

noise_seed

Seed for numpy.random. Used to generate the noise vector if adding noise to the input visibility data. Defaults to 742123.

Type:

float, optional

use_intrinsic_noise_fitting

Fit for the noise level. Defaults to False.

Type:

bool, optional

Instrument Model Parameters
model_instrument

Forward model an instrument (True) or exclude instrumental effects (False). Defaults to True.

Type:

bool, optional

inst_model

Path to a directory containing the instrument model. This directory must at least contain two numpy-compatible files: uvw_model.npy containing the sampled (u, v, w) coordinates with shape (nt, nbls, 3) where nbls is the number of baselines, and redundancy_model.npy containing the number of baselines in a each sampled (u, v, w) with shape (nt, nbls, 1). Please see bayeseor.model.instrument for more details. Used only if include_instrumental_effects is True. Defaults to None.

Type:

str, optional

telescope_latlonalt

Telescope location in latitude (deg), longitude (deg), and altitude (meters). Passed as a list of floats, e.g. ‘[30.1,125.6,80.4]’. Do not put spaces after commas. Required if include_instrumental_effects is True. Defaults to None.

Type:

list of float, optional

drift_scan

Model the instrument in drift scan mode (True) or in phased mode (False). Used only if include_instrumental_effects is True. Defaults to True.

Type:

bool, optional

beam_type

Path to a pyuvdata-compatible beam file or one of ‘uniform’, ‘gaussian’, ‘airy’, ‘gausscosine’, or ‘taperairy’. Used only if include_instrumental_effects is True. Defaults to None.

Type:

str

fwhm_deg

Full width at half maximum of beam in degrees. Used only if include_instrumental_effects is True and beam_type is ‘airy’, ‘gaussian’, or ‘gausscosine’. Defaults to None.

Type:

float, optional

antenna_diameter

Antenna (aperture) diameter in meters. Used only if include_instrumental_effects is True and beam_type is ‘airy’, ‘gaussian’, or ‘gausscosine’. Defaults to None.

Type:

float, optional

cosfreq

Cosine frequency if using a ‘gausscosine’ beam. Used only if include_instrumental_effects is True. Defaults to None.

Type:

float, optional

achromatic_beam

Force the beam to be achromatic. The frequency at which the beam will be calculated is set via beam_ref_freq. Used only if include_instrumental_effects is True. Defaults to False.

Type:

bool, optional

beam_ref_freq

Beam reference frequency in hertz. Used only if include_instrumental_effects is True. Defaults to freq_min.

Type:

bool, optional

beam_peak_amplitude

Peak amplitude of the beam. Used only if include_instrumental_effects is True. Defaults to 1.0.

Type:

float, optional

beam_center

Beam center offsets from the phase center in RA and Dec in degrees. Default behavior is the beam center aligns with the phase center. Passed as a list of floats, e.g. ‘[-1.3,0.01]’. Do not include a space after the comma. Defaults to None.

Type:

list of float, optional

Subharmonic Grid Parameters
nu_sh

Number of pixels on a side for the u-axis in the subharmonic grid model uv-plane. Defaults to None.

Type:

int, optional

nv_sh

Number of pixels on a side for the v-axis in the subharmonic grid model uv-plane. Defaults to nu_sh.

Type:

int, optional

fit_for_shg_amps

Fit for the amplitudes of the subharmonic grid pixels. Defaults to False.

Type:

bool, optional

Tapering Parameters
taper_func

Taper function applied to the frequency axis of the visibilities. Can be any valid argument to scipy.signal.windows.get_window. Defaults to None.

Type:

str, optional

Input Data Parameters
data_path

Path to either a pyuvdata-compatible visibility file or a preprocessed numpy-compatible visibility vector in units of mK sr.

Type:

str

noise_data_path

Path to a preprocessed numpy-compatible noise visibility vector in units of mK sr. Required if calc_noise is False and sigma is None.

Type:

str, optional

ant_str

Antenna downselect string. If data_path points to a pyuvdata-compatible visibility file, ant_str determines what baselines to keep in the data vector. Please see pyuvdata.UVData.select for more details. Defaults to ‘cross’ (cross-correlation baselines only).

Type:

str, optional

bl_cutoff

Baseline length cutoff in meters. If data_path points to a pyuvdata-compatible visibility file, bl_cutoff determines the longest baselines kept in the data vector. Defaults to None (keep all baselines).

Type:

float, optional

form_pI

Form pseudo-Stokes I visibilities. Otherwise, use the polarization specified by pol. Used only if data_path points to a pyuvdata-compatible visibility file. Defaults to True.

Type:

bool, optional

pI_norm

Normalization, N, used in forming pseudo-Stokes I from XX and YY via pI = N * (XX + YY). Used only if data_path points to a pyuvdata-compatible visibility file and form_pI is True. Defaults to 1.0.

Type:

float, optional

pol

Case-insensitive polarization string. Can be one of ‘xx’, ‘yy’, or ‘pI’ for XX, YY, or pseudo-Stokes I polarization, respectively. Used only if data_path points to a pyuvdata-compatible visibility file and form_pI is False. Defaults to ‘xx’.

Type:

str, optional

redundant_avg

Redundantly average the data. Used only if data_path points to a pyuvdata-compatible visibility file. Defaults to False.

Type:

bool, optional

uniform_redundancy

Force the redundancy model to be uniform. Used only if data_path points to a pyuvdata-compatible visibility file. Defaults to False.

Type:

bool, optional

phase_time

The time to which the visibilities will be phased as a Julian date. Used only if drift_scan is False. If drift_scan is False and phase_time is None, phase_time will be automatically set to the central time in the data. Used only if data_path points to a pyuvdata-compatible visibility file. Defaults to None.

Type:

float, optional

calc_noise

Calculate a noise estimate from the visibilities via differencing adjacent times per baseline and frequency. Used only if data_path points to a pyuvdata-compatible visibility file. Defaults to False.

Type:

bool, optional

save_vis

Write visibility vector to disk in out_dir. If calc_noise is True, also save the noise vector. Used only if data_path points to a pyuvdata-compatible visibility file. Defaults to False.

Type:

bool, optional

save_model

Write instrument model (antenna pairs, (u, v, w) sampling, and redundancy model) to disk in out_dir. If phase is True, also save the phasor vector. Used only if data_path points to a pyuvdata-compatible visibility file. Defaults to False.

Type:

bool, optional

parse_args(args_str=None)

Parse arguments from sys.argv or args.

Parameters:

args_str (list of str, optional) – Command line arguments as a list of strings, e.g. ‘[”–config”, “example_config.yaml”]’. If None (default), pulls from sys.argv.

Returns:

args – Namespace of parsed arguments.

Return type:

Namespace

BuildMatrices

class bayeseor.matrices.build.BuildMatrices(nu=None, du_eor=None, nv=None, dv_eor=None, nu_fg=None, du_fg=None, nv_fg=None, dv_fg=None, nf=None, neta=None, deta=None, fit_for_monopole=False, use_shg=False, nu_sh=None, nv_sh=None, nq_sh=None, npl_sh=None, fit_for_shg_amps=False, f_min=None, df=None, nq=0, npl=0, beta=None, sigma=None, nside=None, fov_ra_eor=None, fov_dec_eor=None, fov_ra_fg=None, fov_dec_fg=None, simple_za_filter=True, include_instrumental_effects=True, telescope_latlonalt=(0, 0, 0), nt=None, jd_center=None, dt=None, beam_type=None, beam_center=None, achromatic_beam=False, beam_peak_amplitude=1, fwhm_deg=None, antenna_diameter=None, cosfreq=None, beam_ref_freq=None, drift_scan=True, uvw_array_m=None, bl_red_array=None, phasor=None, effective_noise=None, taper_func=None, array_save_directory='./matrices/', use_sparse_matrices=True, Finv_Fprime=True, verbose=False)

Class for handling matrix construction and arithmetic.

Parameters:
  • nu (int) – Number of pixels on a side for the u axis in the model uv-plane. Defaults to None.

  • du_eor (float) – Fourier mode spacing along the u axis in inverse radians of the EoR model uv-plane. Defaults to None.

  • nv (int) – Number of pixels on a side for the v axis in the model uv-plane. Defaults to None.

  • dv_eor (float) – Fourier mode spacing along the v axis in inverse radians of the EoR model uv-plane. Defaults to None.

  • nu_fg (int) – Number of pixels on a side for the u-axis in the FG model uv-plane. Defaults to None.

  • du_fg (float) – Fourier mode spacing along the u axis in inverse radians of the FG model uv-plane. Defaults to None.

  • nv_fg (int) – Number of pixels on a side for the v-axis in the FG model uv-plane. Defaults to None.

  • dv_fg (float) – Fourier mode spacing along the v axis in inverse radians of the FG model uv-plane. Defaults to None.

  • nf (int) – Number of frequency channels. Defaults to None.

  • neta (int) – Number of Line of Sight (LoS) Fourier modes. Defaults to None.

  • deta (float) – Fourier mode spacing along the eta (line of sight, frequency) axis in inverse Hz. Defaults to None.

  • fit_for_monopole (bool, optional) – Fit for (u, v) = (0, 0) (True) or exclude it from the fit (False). Defaults to False.

  • use_shg (bool, optional) – Use the SubHarmonic Grid (SHG) in the model uv-plane. Defaults to False.

  • nu_sh (int, optional) – Number of pixels on a side for the u-axis in the subharmonic model uv-plane. Defaults to None.

  • nv_sh (int, optional) – Number of pixels on a side for the v-axis in the subharmonic model uv-plane. Defaults to None.

  • nq_sh (int, optional) – Number of large spectral scale modes for each pixel in the subharmonic grid. Defaults to None.

  • npl_sh (int, optional) – Number of power law coefficients used in the large spectral scale model for each pixel in the subharmonic grid. Defaults to None.

  • fit_for_shg_amps (bool, optional) – Fit for the amplitudes of the individual SHG pixels per frequency. Defaults to False.

  • f_min (float) – Minimum frequency in megahertz. Defaults to None.

  • df (float) – Frequency channel width in megahertz. Defaults to None.

  • nq (int, optional) – Number of quadratic modes in the Large Spectral Scale Model (LSSM). Defaults to 0.

  • npl (int, optional) – Number of power law coefficients which replace quadratic modes in the LSSM. Defaults to 0.

  • beta (list of float, optional) – Brightness temperature power law spectral index/indices used in the large spectral scale model. Can be a single spectral index, e.g. [2.63], or multiple spectral indices, e.g. [2.63, 2.82]. Defaults to [2.63, 2.82]. Defaults to None.

  • sigma (float) – Standard deviation of the visibility noise. Defaults to None.

  • nside (int) – HEALPix nside parameter. Defaults to None.

  • fov_ra_eor (float) – Field of view in degrees of the RA axis of the EoR sky model. Defaults to None.

  • fov_dec_eor (float) – Field of view in degrees of the DEC axis of the EoR sky model. Defaults to None.

  • fov_ra_fg (float) – Field of view in degrees of the RA axis of the FG sky model. Defaults to None.

  • fov_dec_fg (float) – Field of view in degrees of the DEC axis of the FG sky model. Defaults to None.

  • simple_za_filter (bool, optional) – Filter pixels in the sky model by zenith angle only (circular FoV) independently along the RA and DEC axes (rectangular FoV). This will be deprecated in a future version. Defaults to True.

  • include_instrumental_effects (bool, optional) – Include instrumental effects like frequency dependent (u, v) sampling and the primary beam. Defaults to True.

  • telescope_latlonalt (tuple, optional) – The latitude, longitude, and altitude of the telescope in degrees, degrees, and meters, respectively. Defaults to (0, 0, 0).

  • nt (int) – Number of times. Defaults to None.

  • jd_center (float) – Central time step of the observation in JD2000 format. Defaults to None.

  • dt (float) – Time cadence of observations in seconds. Defaults to None.

  • beam_type (string) – Beam type to use. Can be ‘uniform’, ‘gaussian’, ‘airy’, ‘taperairy’, or ‘gausscosine’. Defaults to None.

  • beam_center (tuple of floats) – Beam center in (RA, DEC) coordinates and units of degrees. Assumed to be an tuple of offsets along the RA and DEC axes relative to the pointing center of the sky model determined from the instrument model parameters telescope_latlonalt and jd_center. Defaults to None.

  • achromatic_beam (bool, optional) – Force the beam to be achromatic using beam_ref_freq as the reference frequency. Defaults to False.

  • beam_peak_amplitude (float, optional) – Peak amplitude of the beam. Defaults to 1.

  • fwhm_deg (float) – Full Width at Half Maximum (FWHM) of the beam if using a Gaussian beam, or the effective FWHM of the main lobe of an Airy beam from which the diameter of the aperture is calculated. Defaults to None.

  • antenna_diameter (float) – Antenna (aperture) diameter in meters. Used in the calculation of an Airy beam pattern or when using a Gaussian beam with a FWHM that varies as a function of frequency. The FWHM evolves according to the effective FWHM of the main lobe of an Airy beam. Defaults to None.

  • cosfreq (float) – Cosine frequency in radians if using a ‘gausscosine’ beam. Defaults to None.

  • beam_ref_freq (float, optional) – Beam reference frequency in megahertz. Used to fix the beam to be achromatic. Defaults to None.

  • drift_scan (bool, optional) – If True, model a drift scan primary beam, i.e. the beam center drifts across the image space model with time. Defaults to True.

  • uvw_array_m (numpy.ndarray) – Array containing the (u(t), v(t), w(t)) coordinates of the instrument model with shape (nt, nbls, 3). Defaults to None.

  • bl_red_array (numpy.ndarray) – Array containing the number of redundant baselines at each (u(t), v(t), w(t)) in the instrument model with shape (nt, nbls, 1). Defaults to None.

  • phasor (numpy.ndarray, optional) – Array with shape (ndata,) that contains the phasor term used to phase visibilities after performing the nuDFT from HEALPix (l, m, f) to instrumentally sampled, unphased (u, v, f). Defaults to None, i.e. modelling unphased visibilities. Defaults to None.

  • effective_noise (numpy.ndarray) – If the data vector being analyzed contains signal + noise, the effective_noise vector contains the estimate of the noise in the data vector. Must have the shape and ordering of the data vector, i.e. (ndata,). Defaults to None.

  • taper_func (str, optional) – Tapering function to apply to the frequency axis of the model visibilities. Can be any valid argument to scipy.signal.windows.get_window. Defaults to None.

  • array_save_directory (str, optional) – Path to the directory where arrays will be saved. Defaults to ‘./matrices/’.

  • use_sparse_matrices (bool, optional) – Use sparse matrices in place of numpy arrays where possible for memory efficiency. Defaults to True.

  • Finv_Fprime (bool, optional) – Construct the matrix product Finv_Fprime in place from the dense matrices comprising Finv and Fprime to minimize the memory and time required to build the matrix stack (True). In this case, only the matrix product Finv_Fprime is written to disk. Otherwise (False), construct Finv and Fprime independently and save both matrices to disk. Defaults to True.

  • verbose (bool, optional) – Verbose output. Defaults to False.

build_Finv()

Build a multi-frequency NUDFT matrix for image to measurement space.

\(\mathbf{F}^{-1}\) is a a non-uniform DFT matrix that takes a vector of \((l, m, n)\) syk model pixel amplitudes and:

  1. Applies a beam per time and frequency via multi_chan_beam

  2. Transforms to instrumentally sampled, unphased \((u(f), v(f), w(f))\) coordinates from the instrument model

  3. If modelling phased visibilities, applies a phasor vector from the instrument model to phase the visibilities

Notes

  • Finv has shape (ndata, npix_fg * nf).

build_Finv_Fprime()

Build the matrix product \(\mathbf{F}^{-1}\mathbf{F}'\).

This function has been added to reduce the RAM required to build \(\mathbf{F}^{-1}\) and \(\mathbf{F}'\) when the foreground sky model FoV is large. Building the matrix product \(\mathbf{F}^{-1}\mathbf{F}'\) directly from their dense blocks requires much less memory than constructing \(\mathbf{F}^{-1}\) and \(\mathbf{F}'\) independently. If more than one thread is available and threading is available within numpy, this method also significantly reduces the time required to compute the matrix product \(\mathbf{F}^{-1}\mathbf{F}'\).

Notes

  • Finv_Fprime has shape (ndata, (nu*nv - 1)*nf) if modelling the EoR only, (ndata, (nu*nv - 1 + nu_fg*nv_fg - (not fit_for_monopole))*nf) if modelling EoR + foregrounds, (ndata, (nu*nv - 1 + nu_fg*nv_fg - (not fit_for_monopole))*nf + nuv_sh*(nq_sh + fit_for_shg_amps)) if modelling EoR + foregrounds + the subharmonic grid (SHG). Please note that the SHG is not currently supported (see issue #50).

build_Fprime()

Build a multi-frequency NUIDFT matrix for uv to image space.

\(\mathbf{F}'\) takes a rectilinear \((u, v)\) model as a channel-ordered vector and transforms it to HEALPix sky model \((l, m)\) space. Fprime is constructed as a block-diagonal matrix with blocks for the EoR and FG models.

Notes

  • Fprime has shape (npix_fg * nf, (nu*nv - 1) * nf + (nu_fg*nv_fg - (not fit_for_monopole)) * nf.

build_Fprime_Fz()

Build \(\mathbf{F}'\mathbf{F}_z\).

Notes

  • Fprime_Fz has shape (npix_fg*nf, (nu*nv - 1)*(neta - 1) + (nq + 1)*(nu_fg*nv_fg - 1) + fit_for_monopole*(neta + nq)).

build_Fz()

Build a block-diagonal IDFT matrix for eta to frequency space.

\(\mathbf{F}_z\) is constructed as a block-diagonal matrix. Each block is a 1D IDFT matrix which takes a vis-ordered eta space vector and transforms it to a chan-ordered frequency space data vector.

Notes

  • Fz has shape (nf*(nu*nv - 1 + nu_fg*nv_fg - (not fit_for_monopole)), (nu*nv - 1)*(neta - 1) + (nq + 1)*(nu_fg*nv_fg - 1) + fit_for_monopole*(neta + nq)).

build_N()

Build a noise covariance matrix, \(\mathbf{N}\).

Notes

  • N has shape (ndata, ndata).

build_Ninv()

Build an inverse noise covariance matrix, \(\mathbf{N}^{-1}\).

Notes

  • Ninv has shape (ndata, ndata).

build_Ninv_T()

Build \(\mathbf{N}^{-1}\mathbf{T}\).

\(\mathbf{N}^{-1}\mathbf{T}\) computes the inverse noise-covariance-weighted visibility vector from a joint EoR and foreground model vector.

Notes

  • Ninv_T has shape (ndata, (nu*nv - 1)*(neta - 1) + (nq + 1)*(nu_fg*nv_fg - 1) + fit_for_monopole*(neta + nq)).

build_T()

Build \(\mathbf{T}=\mathbf{F}^{-1}\mathbf{F}'\mathbf{F}_z\).

\(\mathbf{T}\) takes a joint EoR and foreground model vector and transforms it to

  1. \((u, v, f)\) space via \(\mathbf{F}_z\)

  2. image space via \(\mathbf{F}'\)

  3. measurement space via \(\mathbf{F}^{-1}\)

Notes

  • T has shape (ndata, (nu*nv - 1)*(neta - 1) + (nq + 1)*(nu_fg*nv_fg - 1) + fit_for_monopole*(neta + nq)).

build_T_Ninv_T()

Build \(\mathbf{T}^\dagger\mathbf{N}^{-1}\mathbf{T}\).

Notes

  • T_Ninv_T is a square matrix with dimension (nu*nv - 1)*(neta - 1) + (nq + 1)*(nu_fg*nv_fg - 1) + fit_for_monopole*(neta + nq).

build_block_T_Ninv_T()

Constructs a block-diagonal representation of T_Ninv_T.

Only used if self.use_instrumental_effects is False.

build_gridding_matrix_co2vo()

Build a chan-ordered to vis-ordered gridding matrix.

This matrix is the transposition of gridding_matrix_vo2co.

Notes

  • gridding_matrix_co2vo is a square matrix with dimension (nu*nv - 1) * nf.

build_gridding_matrix_vo2co()

Build a vis-ordered to chan-ordered gridding matrix for the EoR model.

The gridding matrix takes a \((u, v, f)\)-space vector that is vis-ordered:

  • the first neta entries correspond to the spectrum of the zeroth index model \((u, v)\) pixel

  • the second neta entries correspond to the spectrum of the first index model \((u, v)\) pixel

  • etc.

and converts it to chan ordered:

  • the first nuv entries correspond to the values of the model \((u, v)\) plane at the zeroth frequency channel

  • the second nuv entries correspond to the values of the model \((u, v)\) plane at the first frequency channel

  • etc.

Notes

  • Used for the EoR model in Fz.

  • gridding_matrix_vo2co is a square matrix with dimension (nu*nv - 1) * nf.

build_gridding_matrix_vo2co_fg()

Build a vis-ordered to chan-ordered gridding matrix for the FG model.

The gridding matrix takes a \((u, v, f)\)-space vector that is vis-ordered:

  • the first neta entries correspond to the spectrum of the zeroth index model \((u, v)\) pixel

  • the second neta entries correspond to the spectrum of the first index model \((u, v)\) pixel

  • etc.

and converts it to chan ordered:

  • the first nuv entries correspond to the values of the model \((u, v)\) plane at the zeroth frequency channel

  • the second nuv entries correspond to the values of the model \((u, v)\) plane at the first frequency channel

  • etc.

Notes

  • Used for the FG model in Fz.

  • gridding_matrix_vo2co_fg is a square matrix with dimension (nu_fg*nv_fg - (not fit_for_monopole)) * nf.

build_idft_array_1d()

Build an IDFT matrix for eta to frequency space.

Constructs one block within multi_vis_idft_array_1d.

Notes

  • Used for the EoR model in Fz.

  • idft_array has shape (nf, neta - 1).

  • Excludes \(\eta=0\) which belongs to the FG model.

build_idft_array_1d_fg()

Build a block-diagonal IDFT matrix for eta to frequency space.

Notes

  • Used for the FG model in Fz.

  • idft_array_1d_fg has shape ((nu_fg*nv_fg - (not fit_for_monopole))*nf + nf, (nq + 1)*(nu_fg*nv_fg - 1) + fit_for_monopole*(neta + nq)).

build_idft_array_1d_sh()

Build a block-diagonal IDFT matrix for eta to frequency space.

idft_array_1d_sh is constructed as a block-diagonal matrix. Each block is a 1D IDFT matrix for the eta spectrum of each \((u, v)\) pixel in the SubHarmonic Grid (SHG) model uv-plane.

Notes

  • Used for the SHG model in Fz.

build_matrix_if_it_doesnt_already_exist(matrix_name)

Constructs a matrix with name matrix_name if it doesn’t exist.

This function doesn’t return anything. It instead calls the corresponding build matrix function.

Parameters:

matrix_name (str) – Name of matrix.

build_minimum_sufficient_matrix_stack(clobber_matrices=False, force_clobber=False)

Construct a minimum sufficient matrix stack needed to run BayesEoR.

Parameters:
  • clobber_matrices (bool) – If True, overwrite the existing matrix stack.

  • force_clobber (bool) – If True, delete the old matrix stack without user input. Otherwise, prompt the user to specify wether the matrix stack should be deleted (‘y’) or archived (‘n’).

build_multi_chan_beam()

Build a matrix contating image space beam amplitudes.

Each block-diagonal entry contains the beam amplitude at each HEALPix sampled \((l(t), m(t), n(t))\) for a single time and frequency. Each stack contains nf block-diagonal entries containing the beam amplitudes at all frequencies for a single time.

Notes

  • Used to construct Finv.

  • multi_chan_beam has shape (npix_fg * nf * nt, npix_fg * nf).

build_multi_chan_nudft()

Build a multi-frequency NUDFT matrix for image to measurement space.

Each block in this block-diagonal matrix transforms a set of time-dependent image-space \((l(t), m(t), n(t))\) HEALPix coordinates to unphased, instrumentall sampled, frequency dependent \((u(f), v(f), w(f))\).

Notes

  • Used to construct Finv.

  • multi_chan_nudft has shape (ndata, npix_fg * nf * nt).

  • If use_nvis_nt_nchan_ordering is True: model visibilities will be ordered (nvis*nt) per chan for all channels (old default).

  • If use_nvis_nchan_nt_ordering is True: model visibilities will be ordered (nvis*nchan) per time step for all time steps. This ordering is required when using a drift scan primary beam (current default).

build_multi_chan_nuidft()

Build a multi-frequency NUIDFT matrix for uv to image space.

multi_chan_nuidft is constructed as a block-diagonal matrix. Each block is constructed via build_nuidft_array and represents a 2D non-uniform DFT matrix from rectilinear \((u, v)\) to HEALPix \((l, m)\) for a single frequency.

Notes

  • Used for the EoR model in Fprime.

  • multi_chan_nuidft has shape (npix_fg * nf, (nu*nv - 1) * nf).

build_multi_chan_nuidft_fg()

Build a multi-frequency NUIDFT matrix for uv to image space.

multi_chan_nuidft_fg is constructed as a block-diagonal matrix. Each block is a 2D non-uniform DFT matrix from rectilinear \((u, v)\) to HEALPix \((l, m)\) at a single frequency.

Notes

  • Used for the FG model in Fprime.

  • multi_chan_nuidft_fg has shape (npix_fg * nf, (nu_fg*nv_fg - (not fit_for_monopole)) * nf).

build_multi_vis_idft_array_1d()

Build a block-diagonal IDFT matrix from eta to frequency space.

Notes

  • Used for the EoR model in Fz.

  • multi_vis_idft_array_1d has shape ((nu*nv - 1) * nf, (nu*nv - 1) * (neta - 1)).

build_nuidft_array()

Build a NUIDFT matrix for uv to image space.

This matrix forms a block in multi_chan_nuidft and transforms the EoR model uv-plane to image space at a single frequency. Specifically, nuidft_array transforms a rectilinear \((u, v)\) grid to HEALPix sampled \((l, m)\). The model uv-plane has w=0, so no w or n terms are included in this transformation.

Notes

  • Used for the EoR model in Fprime.

  • nuidft_array has shape (npix_fg, nu*nv - 1).

  • If the EoR and FG models have different FoV values, nuidft_array is reshaped to match the dimensions of the FG model (FoV_FG >= FoV_EoR). The HEALPix pixel ordering must be preserved in this reshaping so that shared pixels between the EoR and FG models are summed together in image-space.

build_nuidft_array_sh()

Build a multi-frequency NUIDFT matrix for uv to image space.

nuidft_array_sh is constructed as a block diagonal matrix. Each block transforms the SubHarmonic Grid (SHG) model uv-plane to HEALPix sampled \((l, m)\) at a single frequency.

Notes

  • Used for the SHG model in Fprime.

  • nuidft_array_sh has shape (npix*nf, nuv_sh*fit_for_shg_amps + nuv_sh*nq_sh).

build_phasor_matrix()

Build a diagonal matrix used to phase visibilities.

The phasor matrix is multiplied elementwise into the visibility vector from \(\mathbf{F}^{-1}\), constructed using unphased \((u, v, w)\) coordinates, to produce phased visibilities.

The phasor matrix is constructed as a diagonal matrix of \(\exp\left[i\,\phi(u(t),v(t),w(t))\right]\) phasor terms from the optional phasor vector in the instrument model.

Notes

  • Used to construct Finv if modelling phased visibilities.

  • phasor_matrix has shape (ndata, ndata).

  • This function assumes that use_nvis_nchan_nt_ordering = True.

build_taper_matrix()

Build a diagonal matrix containing a frequency taper function.

The taper matrix is a diagonal matrix containing a tapering function on the diagonal multiplied elementwise into the visibility vector.

Notes

  • Used to construct Finv if using a taper function.

  • taper_matrix has shape (ndata, ndata).

  • This function assumes that use_nvis_nchan_nt_ordering = True.

check_for_prerequisites(parent_matrix)

Check if parent_matrix requires any child matrices to be built.

Parameters:

parent_matrix (str) – Name of matrix.

Returns:

prerequisites_status – Dictionary containing any required child matrices as keys and integers as values specifying whether each child matrix has been constructed (1 for dense, 2 for sparse) or not (0). If no child matrices are required to build parent_matrix, prerequisites_status is returned as an empty dictionary.

Return type:

dict

check_if_matrix_exists(matrix_name)

Check is hdf5 or npz file with matrix_name exists.

Parameters:

matrix_name (str) – Name of matrix.

Returns:

matrix_available – If matrix exists, 1 or 2 if matrix is an hdf5 or npz file, respectively.

Return type:

int

delete_old_matrix_stack(path_to_old_matrix_stack, confirm_deletion)

Delete or archive an existing matrix stack.

Parameters:
  • path_to_old_matrix_stack (str) – Path to the existing matrix stack.

  • confirm_deletion (str) – If ‘y’, delete existing matrix stack. Otherwise, archive the matrix stack.

dot_product(matrix_a, matrix_b)

Calculate the dot product of matrix_a and matrix_b.

Uses sparse or dense matrix algebra based upon the type of each matrix.

Parameters:
  • matrix_a (array) – First argument.

  • matrix_b (array) – Second argument.

Returns:

ab – Matrix dot product of matrix_a and matrix_b.

Return type:

array

Notes

For dot products of sparse and dense matrices: * dot(sparse, dense) = dense * dot(dense, sparse) = dense * dot(sparse, sparse) = sparse

load_prerequisites(matrix_name)

Load or build any prerequisites for a given matrix.

Parameters:

matrix_name (str) – Name of matrix.

Returns:

prerequisite_matrices_dictionary – Dictionary containing any and all loaded matrix prerequisites.

Return type:

dict

output_data(arr, arr_name, out_dir=None, dataset_name=None)

Write matrix to disk.

Checks if the data is a dense or sparse matrix and calls the corresponding output method.

Parameters:
  • arr (array_like or scipy.sparse) – Array to be written to disk.

  • arr_name (str) – Filename for arr.

  • out_dir (pathlib.Path or str, optional) – Directory in which to write arr. Defaults to self.array_dir.

  • dataset_name (str, optional) – If saving as hdf5, the key used to access arr. Defaults to arr_name.

prepare_matrix_stack_for_deletion(src, clobber_matrices)

Archive an existing matrix stack on disk.

Prepends 'delete_' to the child directory of the matrix stack to be archived.

Parameters:
  • src (str) – Path to existing matrix stack directory.

  • clobber_matrices (bool) – If True, overwrite a previously archived matrix stack.

Returns:

dst – If clobber_matrices is True, path to matrix stack directory to be deleted.

Return type:

str

read_data(matrix_name)

Read matrix from disk as dense or sparse matrix.

This function searches self.array_dir for a file with a prefix matching matrix_name. It first searches for a sparse matrix, i.e. matrix_name + '.npz', and then searches for a dense matrix if a sparse matrix is not found, i.e. matrix_name + '.h5'. If no valid file is found, an error is raised.

Parameters:

matrix_name (str) – Matrix name. For example, to load Finv from disk, matrix_name would be ‘Finv’.

Returns:

data – Loaded matrix as a dense or sparse array.

Return type:

numpy.ndarray or scipy.sparse

read_data_from_hdf5(file_path, dataset_name)

Read array from HDF5 file.

Parameters:
  • file_path (str) – Path to array file.

  • dataset_name (str) – If reading an hdf5 file, the key used to access the dataset.

read_data_from_npz(file_path)

Read sparse matrix from npz file.

Parameters:

file_path (str) – Path to array file.

Notes

  • To maintain sparse matrix attributes, use sparse.load_npz rather than numpy.loadz.

read_data_s2d(matrix_name)

Read matrix from disk and, if necessary, convert to dense matrix.

Parameters:

matrix_name (str) – Matrix name. For example, to load Finv from disk, matrix_name would be ‘Finv’.

sd_block_diag(block_matrices_list)

Generate a block diagonal matrix from a list of matrices.

Parameters:

block_matrices_list (list) – List of input matrices.

Returns:

block_diag_matrix – If self.use_sparse_matrices is True, return a sparse matrix. Otherwise, return a dense numpy array.

Return type:

array

sd_diags(diagonal_vals)

Generate a diagonal matrix from a list of entries in diagonal_vals.

Parameters:

diagonal_vals (array) – Input values to be placed on the matrix diagonal.

Returns:

diagonal_matrix – If self.use_sparse_matrices is True, return a sparse matrix. Otherwise, return a dense numpy array.

Return type:

array

sd_hstack(matrices_list)

Generate a horizontally stacked matrix from a list of matrices.

Parameters:

matrices_list (list) – List of input matrices.

Returns:

hstack_matrix – If self.use_sparse_matrices is True, return a sparse matrix. Otherwise, return a dense numpy array.

Return type:

array

sd_vstack(matrices_list)

Generate a vertically stacked matrix from a list of matrices.

Parameters:

matrices_list (list) – List of input matrices.

Returns:

vstack_matrix – If self.use_sparse_matrices is True, return a sparse matrix. Otherwise, return a dense numpy array.

Return type:

array

sparse2ndarray(arr)

Convert scipy.sparse matrix to dense numpy.ndarray.

Parameters:

arr (scipy.sparse) – Sparse matrix.

Returns:

arr_ndarray – Dense representation of arr.

Return type:

numpy.ndarray

write_version_info()

Write version info to disk.

Cosmology

class bayeseor.cosmology.Cosmology

Class for performing cosmological distance calculations using astropy.cosmology.Planck18.

dL_df(z)

Comoving differential distance at redshift per frequency.

Parameters:

z (float) – Input redshift.

Returns:

dl_df – Conversion factor relating a bandwidth in Hz to a comoving size in Mpc at redshift z.

Return type:

float

dL_dth(z)

Comoving transverse distance per radian in Mpc.

Parameters:

z (float) – Input redshift.

Returns:

dl_dth – Conversion factor relating an angular size in radians to a comoving transverse size in Mpc at redshift z.

Return type:

float

f2z(f)

Convert a frequency f in Hz to redshift relative to self.f_21.

Parameters:

f (float) – Input frequency in Hz.

Returns:

z – Redshift corresponding to frequency f.

Return type:

float

inst_to_cosmo_vol(z)

Conversion factor to go from an instrumentally sampled volume in sr Hz to a comoving cosmological volume in Mpc^3.

Parameters:

z (float) – Input redshift.

Returns:

i2cV – Volume conversion factor for sr Hz –> Mpc^3 at redshift z.

Return type:

float

z2f(z)

Convert a redshift z relative to self.f_21 to a frequency in Hz.

Parameters:

z (float) – Input redshift.

Returns:

f – Frequency corresponding to redshift z.

Return type:

float

DataContainer

class bayeseor.analyze.analyze.DataContainer(dirnames, dir_prefix=None, sampler='multinest', store_samples=False, cred_intervals=[68, 95], calc_uplims=True, uplim_quantile=0.95, uplim_inds=None, posterior_weighted=True, Nhistbins=31, density=False, calc_kurtosis=False, ps_kind='dmps', temp_unit='mK', little_h_units=False, expected_ps=None, expected_dmps=None, labels=None)

Class for reading and analyzing files output by BayesEoR.

Parameters:
  • dirnames (array-like of str) – Array-like of BayesEoR output directory names.

  • dir_prefix (str, optional) – Prefix to append to dirnames. Defaults to None, i.e. the entries in dirnames are assumed to be valid paths to BayesEoR output directories.

  • sampler (str, optional) – Case insensitive sampler name, e.g. ‘MultiNest’ or ‘multinest’ (default). The only currently supported sampler is MultiNest (default).

  • store_samples (bool, optional) – If store_samples is True, store all samples as an attribute, self.samples. Defaults to False.

  • cred_intervals (float or list of float, optional) – Credibility intervals as percentages used for uncertainty calculations and errorbars on plots. Defaults to [68, 95], i.e. compute 68% and 95% credibility intervals.

  • calc_uplims (bool, optional) – If calc_uplims is True, calculate the upper limit of each k bin’s posterior as the 95th percentile. Defaults to False.

  • uplim_quantile (float, optional) – Quantile in [0, 1]. Defaults to 0.95. Only used if calc_uplims is True.

  • uplim_inds (array-like, optional) – Array-like of True for non-detections and False for detections. Can have shape (Nkbins,), where Nkbins is the number of spherically averaged k bins with power spectrum posteriors in the sampler output or shape (len(dirnames), Nkbins). If Nkbins varies in each file, each entry in uplims must have len(uplims[i]) == Nkbins for that particular file.

  • posterior_weighted (bool, optional) – If posterior_weighted is True, use the joint posterior probability as weights in the calculation of the individual power spectrum coefficient posteriors and the associated estimates, uncertainties, and upper limits. Defaults to False.

  • Nhistbins (int, optional) – Number of histogram bins for each k bin’s posterior distribution.

  • density (bool, optional) – If density is True, compute the posterior as a probability density function. Defaults to False, i.e. plot counts.

  • calc_kurtosis (bool, optional) – If calc_kurtosis is True, calculate the kurtosis of each k bin’s posterior distribution.

  • ps_kind (str, optional) – Case insensitive power spectrum kind in the sampler output file. Can be ‘ps’ or ‘dmps’ for the power spectrum, \(P(k)\), or dimensionless power spectrum, \(\Delta^2(k)\). Defaults to ‘dmps’.

  • temp_unit (str, optional) – Either ‘mK’ or ‘K’. The temperature unit of the power spectrum. The default output from BayesEoR is ‘mK’ (default).

  • little_h_units (bool, optional) – If little_h_units is True, power spectrum samples are assumed to be in units of little h. Defaults to False. Currently, BayesEoR picks the Planck 2018 value of the Hubble constant but we plan to output power spectra in little h units in the future. This kwarg has thus been added for future use.

  • expected_ps (float or array-like, optional) – Expected power spectrum, \(P(k)\). Can be a single float (for a flat \(P(k)\)) or an array-like with length equal to the number of spherically averaged k bins.

  • expected_dmps (float or array-like, optional) – Expected dimensionless power spectrum, \(\Delta^2(k)\). Can be a single float or an array-like with length equal to the number of spherically averaged k bins.

  • labels (array-like of str, optional) – Array-like containing strings with shorthand labels for each directory in dirnames. Used in figure legends for plotting functions. Defaults to None (no labels).

calculate_expected_ps(expected_ps=None, expected_dmps=None)

Calculated the expected power spectrum for each output file.

Sets either self.expected_ps if self.ps_kind is ‘ps’ or self.expected_dmps if self.ps_kind is ‘dmps’.

Parameters:
  • expected_ps (float or array-like, optional) – Expected power spectrum, \(P(k)\). Can be a single float (for a flat \(P(k)\)) or an array-like with length equal to the number of spherically averaged k bins.

  • expected_dmps (float or array-like, optional) – Expected dimensionless power spectrum, \(\Delta^2(k)\). Can be a single float or an array-like with length equal to the number of spherically averaged k bins.

dmps_to_ps(dmps, ks)

Convert \(\Delta^2(k)\) to \(P(k)\).

Parameters:
  • dmps (float or numpy.ndarray) – Float or ndarray of floats containing the dimensionless power spectrum amplitude(s). If an ndarray, must have the same shape as ks.

  • ks (numpy.ndarray) – Array of k values.

get_posterior_data(fp, Nkbins, posterior_weighted=False, cred_intervals=[68, 95], uplim_quantile=0.95, Nhistbins=31, density=False, log_priors=True, return_samples=False)

Load sampler output and form posteriors for each k bin.

Parameters:
  • fp (str or Path) – Path to sampler output file.

  • Nkbins (int) – Number of spherically averaged k bins.

  • posterior_weighted (bool, optional) – If posterior_weighted is True, use the joint posterior probability as weights in the calculation of the individual power spectrum coefficient posteriors and the associated estimates, uncertainties, and upper limits. Defaults to False.

  • cred_intervals (float or list of floats, optional) – Credibility intervals as percentages. Defaults to [68, 95], i.e. compute 68% and 95% credibility intervals.

  • uplim_quantile (float, optional) – Quantile in [0, 1]. Defaults to 0.95. Only used if self.calc_uplims is True.

  • Nhistbins (int, optional) – Number of histogram bins for each k bin’s posterior distribution. Defaults to 31.

  • density (bool, optional) – If density is True, compute the posterior as a probability density function. Defaults to False, i.e. plot counts.

  • log_priors (bool, optional) – If log_priors is True (default), power spectrum samples are assumed to be in log10 units and are first linearized prior to calculating e.g. upper limits. Otherwise, the power spectrum samples are assumed to be in linear units.

  • return_samples (bool, optional) – If return_samples is True, return all samples. Defaults to False.

Returns:

  • posteriors (numpy.ndarray) – Posterior distributions for each k bin with shape (Nkbins, Nhistbins) where Nkbins is the number of spherically averaged k bins in fp.

  • avgs (numpy.ndarray) – Average of each k bin.

  • medians (numpy.ndarray) – Median value of each k bin.

  • ci_dict (dict) – Dictionary with the credibility interval(s) as key(s) and nested dictionaries for each credibility interval indexed by ‘lo’ and ‘hi’ for the low and high bounds of the credibility interval, respectively.

  • uplims (numpy.ndarray) – uplim_quantile-th quantile of each k bin. Returned only if self.calc_uplims is True.

  • kurtoses (numpy.ndarray) – Kurtosis of each k bin’s posterior. Returned only if self.calc_kurtosis is True.

  • samples (numpy.ndarray) – Samples for each power power spectrum coefficient with shape (Nsamples, Nkbins + 2). Returned only if return_samples is True.

plot_posteriors(plot_height=1.0, plot_width=7.0, hspace=0, colors=None, cmap=None, lw=3, ls_expected='--', log_y=False, ymin=1e-16, show_k_vals=True, plot_priors=False, suptitle=None, fig=None, axs=None)

Plot power spectrum posteriors.

Parameters:
  • plot_height (float, optional) – Subplot height. Defaults to 1.0.

  • plot_width (float, optional) – Subplot width. Defaults to 7.0.

  • hspace (float, optional) – matplotlib gridspec height space. Defaults to 0.05. Only used if plot_diff or plot_fracdiff is True.

  • colors (array-like of str, optional) – Array-like of valid matplotlib color strings. Must have length self.Ndirs. If None (default), use the default matplotlib color sequence.

  • cmap (callable matplotlib colormap instance, optional) – Callable matplotlib colormap instance, e.g. matplotlib.pyplot.cm.viridis. If None (default), use the default matplotlib color sequence.

  • lw (float, optional) – matplotlib line width. Defaults to 3.

  • ls_expected (str, optional) – Any valid matplotlib line style string. Used for plotting the expected power spectrum (used only if self.expected_ps or self.expected_dmps is not None). Defaults to ‘–‘.

  • log_y (bool, optional) – If log_y is True, plot the amplitudes of the posteriors in log10 units. Otherwise, plot the amplitudes of the posteriors in linear units (default).

  • ymin (float, optional) – Minimum value for the y axis. Used if log_y is True to avoid plotting very small y values. Defaults to 1e-16.

  • show_k_vals (bool, optional) – If show_k_vals is True (default), print the value of the k bin associated with each posterior in the upper left corner of each posterior subplot.

  • plot_priors (bool, optional) – If plot_priors is True, plot the prior bounds as shaded regions for each k bin. Defaults to False.

  • suptitle (str, optional) – Figure suptitle string. Defaults to None.

  • fig (Figure, optional) – matplotlib Figure instance. Used internally when called by self.plot_power_spectra_and_posteriors.

  • axs (Axes) – matplotlib Axes instance(s). Used internally when called by self.plot_power_spectra_and_posteriors.

plot_power_spectra(cred_interval=68, uplim_inds=None, plot_height=4.0, plot_width=7.0, hspace=0.05, height_ratios=[1, 0.5], x_offset=0, zorder_offset=0, labels=None, colors=None, cmap=None, marker='o', capsize=3, lw=3, ls_expected='--', plot_diff=False, plot_fracdiff=False, ylim=None, ylim_diff=[-1, 1], plot_priors=False, legend_loc='best', legend_ncols=1, figlegend=False, top=0.85, suptitle=None, fig=None, axs=None)

Plot the power spectrum as the median with a credibility interval.

Parameters:
  • cred_interval (float, optional) – Credibility interval as a percentage to plot as the uncertainty. Defaults to 68.

  • uplim_inds (array-like, optional) – Array-like of True for non-detections and False for detections. Can have shape (Nkbins,), where Nkbins is the number of spherically averaged k bins with power spectrum posteriors in the sampler output or shape (len(dirnames), Nkbins). If Nkbins varies in each file, each entry in uplims must have len(uplims[i]) == Nkbins for that particular file. If uplim_inds is None (default), use self.uplim_inds. Otherwise, use uplim_inds in place of self.uplim_inds.

  • plot_height (float, optional) – Subplot height. Defaults to 4.0.

  • plot_width (float, optional) – Subplot width. Defaults to 7.0.

  • hspace (float, optional) – matplotlib gridspec height space. Defaults to 0.05. Only used if plot_diff or plot_fracdiff is True.

  • height_ratios (array-like, optional) – matplotlib gridspec subplot height ratios. Defaults to [1, 0.5], i.e. the top plot will be twice as tall as the bottom subplot. Only used if plot_diff or plot_fracdiff is True.

  • x_offset (float, optional) – x-axis offset for plotting multiple results on a single subplot. If x_offset > 0, data points for different analyses will be offset along the x-axis to better distinguish overlapping data.

  • zorder_offset (int, optional) – matplotlib zorder offset for plotting multiple results on a single subplot. If zorder_offset > 0, data points for different analyses will be offset along the “z” axis (plot data points over or under each other).

  • labels (array-like of str, optional) – Array-like of label strings for each analysis result. If no labels are provided, checks for labels in self.labels. If self.labels is not None and labels is not None, the labels in labels will be used instead of self.labels. Must have length self.Ndirs. Defaults to None, i.e. use self.labels if not None otherwise use no labels.

  • colors (array-like of str, optional) – Array-like of valid matplotlib color strings. Must have length self.Ndirs. If None (default), use the default matplotlib color sequence.

  • cmap (callable matplotlib colormap instance, optional) – Callable matplotlib colormap instance, e.g. matplotlib.pyplot.cm.viridis. If None (default), use the default matplotlib color sequence.

  • marker (str, optional) – matplotlib marker string. Defaults to ‘o’.

  • capsize (float, optional) – Errorbar cap size. Defaults to 3.

  • lw (float, optional) – matplotlib line width. Defaults to 3.

  • ls_expected (str, optional) – Any valid matplotlib line style string. Used for plotting the expected power spectrum (used only if self.expected_ps or self.expected_dmps is not None). Defaults to ‘–‘.

  • plot_diff (bool, optional) – If plot_diff is True and self.expected_ps or self.expected_dmps is not None, plot the difference between each analysis’ power spectrum and the expected power spectrum. If both plot_diff and plot_fracdiff are True, plot_diff will be set to False and the fractional difference will be plotted instead. Defaults to False.

  • plot_fracdiff (bool, optional) – If plot_fracdiff is True and self.expected_ps or self.expected_dmps is not None, plot the fractional difference between each analysis’ power spectrum and the expected power spectrum. If both plot_diff and plot_fracdiff are True, plot_diff will be set to False and the fractional difference will be plotted instead. Defaults to False.

  • ylim (array-like, optional) – matplotlib ylim for the power spectrum subplot. Defaults to None (scales the y axis limits according to the data).

  • ylim_diff (array-like, optional) – matplotlib ylim for the (fractional) difference subplot if plot_diff or plot_fracdiff is True.

  • plot_priors (bool, optional) – If plot_priors is True, plot the prior bounds as shaded regions for each k bin. Defaults to False.

  • legend_loc (str, optional) – Any valid matplotlib legend locator string. Used only if figlegend is False. Defaults to ‘best’.

  • legend_ncols (int, optional) – Number of columns in the legend. Defaults to 1 unless figlegend is True. In the latter case, the default value is set to self.Ndirs + 1.

  • figlegend (bool, optional) – If figlegend is True, use a figure legend instead of an axes legend. Defaults to False.

  • top (float, optional) – Sets the top of the power spectrum subplot in figure fraction units (0, 1]. Defaults to 0.85.

  • suptitle (str, optional) – Figure suptitle string. Defaults to None.

  • fig (Figure, optional) – matplotlib Figure instance. Used internally when called by self.plot_power_spectra_and_posteriors.

  • axs (Axes) – matplotlib Axes instance(s). Used internally when called by self.plot_power_spectra_and_posteriors.

plot_power_spectra_and_posteriors(cred_interval=68, uplim_inds=None, plot_height_ps=4.0, plot_width=7.0, hspace_ps=0.05, height_ratios_ps=[1, 0.5], x_offset=0, zorder_offset=0, labels=None, colors=None, cmap=None, marker='o', capsize=3, lw=3, ls_expected='--', plot_diff=False, plot_fracdiff=False, ylim_ps=None, ylim_diff_ps=[-1, 1], plot_priors=False, plot_height_post=1.0, hspace_post=0.01, log_y=False, ymin_post=1e-16, show_k_vals=True, legend_ncols=0, figlegend=True, top=0.875, right_ps=0.46, left_post=0.54, suptitle=None)

Make a plot containing power spectra and posteriors as subplots.

Calls self.plot_power_spectra and self.plot_posteriors.

Parameters:
  • cred_interval (float, optional) – Credibility interval as a percentage to plot as the uncertainty. Defaults to 68.

  • uplim_inds (array-like, optional) – Array-like of True for non-detections and False for detections. Can have shape (Nkbins,), where Nkbins is the number of spherically averaged k bins with power spectrum posteriors in the sampler output or shape (len(dirnames), Nkbins). If Nkbins varies in each file, each entry in uplims must have len(uplims[i]) == Nkbins for that particular file. If uplim_inds is None (default), use self.uplim_inds. Otherwise, use uplim_inds in place of self.uplim_inds.

  • plot_height_ps (float, optional) – Subplot height for the power spectra subplot(s). Defaults to 4.0.

  • plot_width (float, optional) – Subplot width for both the power spectra subplot(s) and the posterior subplots. Defaults to 7.0.

  • hspace_ps (float, optional) – matplotlib gridspec height space for the power spectra subplot(s). Defaults to 0.05. Only used if plot_diff or plot_fracdiff is True.

  • height_ratios_ps (array-like, optional) – matplotlib gridspec subplot height ratios for the power spectra subplots. Defaults to [1, 0.5], i.e. the top plot will be twice as tall as the bottom subplot. Only used if plot_diff or plot_fracdiff is True.

  • x_offset (float, optional) – x-axis offset for plotting multiple results on a single subplot. If x_offset > 0, data points for different analyses will be offset along the x-axis to better distinguish overlapping data.

  • zorder_offset (int, optional) – matplotlib zorder offset for plotting multiple results on a single subplot. If zorder_offset > 0, data points for different analyses will be offset along the “z” axis (plot data points over or under each other).

  • labels (array-like of str, optional) – Array-like of label strings for each analysis result. If no labels are provided, checks for labels in self.labels. If self.labels is not None and labels is not None, the labels in labels will be used instead of self.labels. Must have length self.Ndirs. Defaults to None, i.e. use self.labels if not None otherwise use no labels.

  • colors (array-like of str, optional) – Array-like of valid matplotlib color strings. Must have length self.Ndirs. If None (default), use the default matplotlib color sequence.

  • cmap (callable matplotlib colormap instance, optional) – Callable matplotlib colormap instance, e.g. matplotlib.pyplot.cm.viridis. If None (default), use the default matplotlib color sequence.

  • marker (str, optional) – matplotlib marker string. Defaults to ‘o’.

  • capsize (float, optional) – Errorbar cap size. Defaults to 3.

  • lw (float, optional) – matplotlib line width. Defaults to 3.

  • ls_expected (str, optional) – Any valid matplotlib line style string. Used for plotting the expected power spectrum (used only if self.expected_ps or self.expected_dmps is not None). Defaults to ‘–‘.

  • plot_diff (bool, optional) – If plot_diff is True and self.expected_ps or self.expected_dmps is not None, plot the difference between each analysis’ power spectrum and the expected power spectrum. If both plot_diff and plot_fracdiff are True, plot_diff will be set to False and the fractional difference will be plotted instead. Defaults to False.

  • plot_fracdiff (bool, optional) – If plot_fracdiff is True and self.expected_ps or self.expected_dmps is not None, plot the fractional difference between each analysis’ power spectrum and the expected power spectrum. If both plot_diff and plot_fracdiff are True, plot_diff will be set to False and the fractional difference will be plotted instead. Defaults to False.

  • ylim_ps (array-like, optional) – matplotlib ylim for the power spectrum subplot. Defaults to None (scales the y axis limits according to the data).

  • ylim_diff_ps (array-like, optional) – matplotlib ylim for the (fractional) difference power spectrum subplot if plot_diff or plot_fracdiff is True.

  • plot_priors (bool, optional) – If plot_priors is True, plot the prior bounds as shaded regions for each k bin. Defaults to False.

  • plot_height_post (float, optional) – Subplot height for the posterior subplots. Defaults to 1.0.

  • hspace_post (float, optional) – matplotlib gridspec height space for the posterior subplots. Defaults to 0.01.

  • log_y (bool, optional) – If log_y is True, plot the amplitudes of the posteriors in log10 units. Otherwise, plot the amplitudes of the posteriors in linear units (default).

  • ymin (float, optional) – Minimum value for the y axis of the posterior subplots. Used if log_y is True to avoid plotting very small y values. Defaults to 1e-16.

  • show_k_vals (bool, optional) – If show_k_vals is True (default), print the value of the k bin associated with each posterior in the upper left corner of each posterior subplot.

  • legend_ncols (int, optional) – Number of columns in the legend. Defaults to self.Ndirs + 1.

  • figlegend (bool, optional) – If figlegend is True (default), use a figure legend. Otherwise, no legend is shown.

  • top (float, optional) – Sets the top of the power spectrum and posterior subplots in figure fraction units (0, 1]. Defaults to 0.875.

  • right_ps (float, optional) – Sets the right edge of the power spectrum subplot(s) in figure fraction units (0, 1]. Defaults to 0.46.

  • left_post (float, optional) – Sets the left edge of the posterior subplot(s) in figure fraction units (0, 1]. Defaults to 0.54.

  • suptitle (str, optional) – Figure suptitle string. Defaults to None.

ps_to_dmps(ps, ks)

Convert \(P(k)\) to \(\Delta^2(k)\).

Parameters:
  • ps (float or numpy.ndarray) – Float or ndarray of floats containing the power spectrum amplitude(s). If an ndarray, must have the same shape as ks.

  • ks (numpy.ndarray) – Array of k values.

GPUInterface

class bayeseor.gpu.GPUInterface(base_dir=None, rank=0, verbose=True)

Class to interface with GPUs.

Parameters:
  • base_dir (str) – Path to directory containing a shared library. Defaults to the parent directory of __file__.

  • rank (int) – MPI rank.

  • verbose (bool) – If True (default), print status.

Healpix

Class containing HEALPix functions for the image domain model.

class bayeseor.model.healpix.Healpix(fov_ra_eor=None, fov_dec_eor=None, fov_ra_fg=None, fov_dec_fg=None, simple_za_filter=False, nside=256, telescope_latlonalt=(-30.72152777777791, 21.428305555555557, 1073.0000000093132), jd_center=None, nt=1, dt=None, beam_type=None, peak_amp=1.0, fwhm_deg=None, diam=None, cosfreq=None, tanh_freq=None, tanh_sl_red=None, pol='xx', freq_interp_kind='cubic')

Class to store and manipulate HEALPix maps using astropy_healpix functions.

Parameters:
  • fov_ra_eor (float) – Field of view in degrees of the RA axis of the EoR sky model.

  • fov_dec_eor (float, optional) – Field of view in degrees of the DEC axis of the EoR sky model. Defaults to fov_ra_eor.

  • fov_ra_fg (float, optional) – Field of view in degrees of the RA axis of the FG sky model.

  • fov_dec_fg (float, optional) – Field of view in degrees of the DEC axis of the FG sky model. Defaults to fov_ra_fg.

  • simple_za_filter (boolean, optional) – If True, filter pixels in the FoV by zenith angle only. Otherwise, filter pixels in a rectangular region set by the FoV values along RA and DEC.

  • nside (int) – Nside resolution of the HEALPix map. Defaults to 256.

  • telescope_latlonalt (tuple, optional) – Tuple containing the latitude, longitude, and altitude of the telescope in degrees, degrees, and meters, respectively. Defaults to the location of the HERA telescope, i.e. (-30.72152777777791, 21.428305555555557, 1073.0000000093132).

  • jd_center (float) – Central time step of the observation in JD2000 format.

  • nt (int, optional) – Number of time integrations. Defaults to 1.

  • dt (float, optional) – Integration time in seconds. Required if nt > 1.

  • beam_type (str, optional) – Can be either a path to a pyuvdata-compatible beam file or one of ‘uniform’, ‘gaussian’, or ‘airy’. Defaults to ‘uniform’.

  • peak_amp (float, optional) – Peak amplitude of the beam. Defaults to 1.0.

  • fwhm_deg (float, optional) – Full width at half maximum of the beam in degrees. Required if beam_type is ‘gaussian’ or ‘gausscosine’.

  • diam (float, optional) – Antenna (aperture) diameter in meters. Used if beam_type is ‘airy’.

  • cosfreq (float, optional) – Cosine frequency in inverse radians. Used if beam_type is ‘gausscosine’.

  • tanh_freq (float, optional) – Exponential frequency (rate parameter) in inverse radians. Used if beam_type is ‘tanhairy’.

  • tanh_sl_red (float, optional) – Airy sidelobe amplitude reduction as a fractional percent. For example, passing 0.99 reduces the sidelobes by 0.01, i.e. two orders of magnitude. Used if beam_type is ‘tanhairy’.

  • pol (str, optional) – Polarization string. Can be ‘xx’, ‘yy’, or ‘pI’. Only used if beam_type is a path to a pyuvdata-compatible beam file. Defaults to ‘xx’.

  • freq_interp_kind (str, optional) – Frequency interpolation kind. Please see scipy.interpolate.interp1d for valid options and more details. Defaults to ‘cubic’.

calc_lmn_from_radec(time, ra, dec, return_azza=False, radec_offset=None)

Return arrays of (l, m, n) coordinates in radians for all (RA, DEC).

Parameters:
  • time (float) – Julian date used in ICRS to AltAz coordinate frame conversion.

  • ra (numpy.ndarray) – RA values in degrees.

  • dec (numpy.ndarray) – DEC values in degrees.

  • return_azza (boolean) – If True, return both (l, m, n) and (az, za) coordinate arrays. Otherwise return only (l, m, n). Defaults to ‘False’.

  • radec_offset (tuple of floats) – Will likely be deprecated.

Returns:

  • l (numpy.ndarray of floats) – Array containing the EW direction cosine of each HEALPix pixel.

  • m (numpy.ndarray of floats) – Array containing the NS direction cosine of each HEALPix pixel.

  • n (numpy.ndarray of floats) – Array containing the radial direction cosine of each HEALPix pixel.

Notes

get_beam_vals(az, za, freq=None)

Get an array of beam values from (az, za) coordinates.

If self.beam_type is ‘gaussian’, this function assumes that the beam width is symmetric along the l and m axes.

Parameters:
  • az (numpy.ndarray of floats) – Azimuthal angle of each pixel in radians.

  • za (numpy.ndarray of floats) – Zenith angle of each pixel in radians.

  • freq (float, optional) – Frequency in Hz.

Returns:

beam_vals – Array containing beam amplitude values at each (az, za).

Return type:

numpy.ndarray

get_extent_ra_dec(fov_ra, fov_dec, fov_fac=1.0)

Get the sampled extent of the sky in RA and DEC.

Parameters:
  • fov_ra (float) – Field of view in degrees of the RA axis.

  • fov_dec (float) – Field of view in degrees of the DEC axis.

  • fov_fac (float) – Scaling factor for the sampled extent.

Returns:

  • range_ra (tuple) – fov_fac scaled (min, max) sampled RA values.

  • range_dec (tuple) – fov_fac scaled (min, max) sampled DEC values.

get_pixel_filter(fov_ra, fov_dec, return_radec=False, inverse=False, simple_za_filter=True)

Return HEALPix pixel indices lying inside an observed region.

This function gets the HEALPix pixel indices for all pixel centers lying inside

  • a rectangle with equal arc length on all sides if simple_za_filter is True

  • a circle with radius fov_ra if simple_za_filter is False

Parameters:
  • fov_ra (float) – Field of view in degrees of the RA axis.

  • fov_dec (float) – Field of view in degrees of the DEC axis.

  • return_radec (bool, optional) – Return the (RA, DEC) coordinates associated with each pixel center. Defaults to False.

  • inverse (boolean, optional) – If False, return the pixels within the observed region. If True, return the pixels outside the observed region.

  • simple_za_filter (boolean, optional) – If True (default), return the pixels inside a circular region defined by za <= fov_ra/2 where za is the zenith angle. Otherwise, return the pixels inside a rectangular region with equal arc length on all sides.

Returns:

  • pix (numpy.ndarray) – HEALPix pixel numbers lying within the observed region set by fov_ra and fov_dec.

  • ra (numpy.ndarray) – RA values for each pixel center. Only returned if return_radec is True.

  • dec (numpy.ndarray) – DEC values for each pixel center. Only returned if return_radec is True.

Notes

  • The rectangular pixel selection (simple_za_filter is False) has been left for posterity. It has been found to be flawed when the field of view along right ascension becomes large (see issue #11 in the BayesEoR repo for more details). We advise setting simple_za_filter to True (the default) to avoid any potential issues with the rectangular pixel selections.

MaximumAPosteriori

class bayeseor.analyze.maxap.MaximumAPosteriori(*, config: Path | str, data_path: str | None = None, array_dir: str | None = None, verbose: bool = False)

Class for performing maximum a posteriori calculations.

Parameters:
  • config (str) – Path to a BayesEoR yaml configuration file.

  • data_path (str, optional) – Path to either a pyuvdata-compatible file containing visibilities or a numpy-compatible file containing a preprocessed visibility vector. Supersedes the data_path in config. Defaults to None (use the data_path in config).

  • array_dir (str, optional) – Path to a directory containing T, Ninv, and T_Ninv_T. Supersedes the array_dir created based on the parameters in config. Defaults to None (use the array_dir created from the parameters in config).

  • verbose (bool) – Verbose output. Defaults to False.

calculate_dmps(ps: float | _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None = None, dmps: float | _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None = None)

Calculated the expected dimensionless power spectrum.

Parameters:
  • ps (float or array-like, optional) – Expected power spectrum, \(P(k)\). Can be a single float (for a flat power spectrum) or an array-like with shape (self.k_vals.size,). Required if dmps is None.

  • dmps (float or array-like, optional) – Expected dimensionless power spectrum, \(\Delta^2(k)\). Can be a single float (for a flat dimensionless power spectrum) or an array-like with shape (self.k_vals.size,). Required if ps is None.

Returns:

dmps – Array of dimensionless power spectrum amplitudes with shape (self.k_vals.size,).

Return type:

numpy.ndarray

calculate_prior_covariance(ps: float | _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None = None, dmps: float | _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None = None)

Calculate the prior covariance matrix \(\Phi^{-1}\).

Parameters:
  • ps (float or array-like, optional) – Expected power spectrum, \(P(k)\). Can be a single float (for a flat power spectrum) or an array-like with shape (self.k_vals.size,). Required if dmps is None.

  • dmps (float or array-like, optional) – Expected dimensionless power spectrum, \(\Delta^2(k)\). Can be a single float (for a flat dimensionless power spectrum) or an array-like with shape (self.k_vals.size,). Required if ps is None.

Returns:

PhiI – Diagonal of the prior covariance matrix.

Return type:

ndarray

extract_bl_vis(vis_vec: ndarray)

Form a dictionary of visibility waterfalls indexed by baseline.

The input visibility vector should have shape (nbls*nf*nt,) where nbls is the number of baselines, nf the number of frequencies, and nt the number of times. The ordering of the visibility vector is expected to follow the ordering of the visibilities in Finv, i.e. the first nbls visibilities correspond to all baselines at the 0th frequency and 0th time, the second nbls visibilities correspond to all baselines at the 1st frequency and 0th time, etc. The baseline axis evolves most rapidly, then frequency, and then time. In this ordering, a visibility waterfall for the i-th indexed baseline, i_bl, with shape (nt, nf) can be obtained via

bl_vis = vis_vec[i_bl::nbls].reshape(nt, nf)
Parameters:

vis_vec (numpy.ndarray) – Visibility vector with shape (nbls*nf*nt,).

Returns:

bl_vis – Dictionary with keys of either baseline index or antenna pair tuple if self.antpairs exists.

Return type:

dict

map_estimate(ps: float | _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None = None, dmps: float | _Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str] | None = None, return_prior_cov: bool = False)

Calculate the maximum a posteriori (MAP) model coefficients.

Must pass one of ps or dmps which is required to calculate the prior covariance matrix.

Parameters:
  • ps (float or array-like, optional) – Expected power spectrum, \(P(k)\). Can be a single float (for a flat power spectrum) or an array-like with shape (self.k_vals.size,). Required if dmps is None.

  • dmps (float or array-like, optional) – Expected dimensionless power spectrum, \(\Delta^2(k)\). Can be a single float (for a flat dimensionless power spectrum) or an array-like with shape (self.k_vals.size,). Required if ps is None.

  • return_prior_cov (bool, optional) – Return the diagonal of the prior covariance matrix. Defaults to False.

Returns:

  • map_coeffs (numpy.ndarray) – MAP coefficients of the model.

  • map_vis (numpy.ndarray) – MAP model visibilities.

  • log_post (float) – Log posterior probability.

  • prior_cov (numpy.ndarray) – Diagonal of the prior covariance matrix. Only returned if return_prior_cov is True.

Notes

  • This function does not currently support noise fitting

PowerSpectrumPosteriorProbability

class bayeseor.posterior.PowerSpectrumPosteriorProbability(T_Ninv_T, dbar, k_vals, k_cube_voxels_in_bin, nuv, neta, nf, nq, Ninv, d_Ninv_d, redshift, ps_box_size_ra_Mpc, ps_box_size_dec_Mpc, ps_box_size_para_Mpc, include_instrumental_effects=True, log_priors=False, uprior_inds=None, use_LWM_Gaussian_prior=False, inverse_LW_power=1e-16, dimensionless_PS=True, block_T_Ninv_T=None, intrinsic_noise_fitting=False, fit_for_spectral_model_parameters=False, pl_max=None, pl_grid_spacing=None, use_shg=False, return_Sigma=False, rank=0, use_gpu=True, verbose=False, print_rate=100, debug=False, print_debug=False)

Class containing posterior calculation functions.

Parameters:
  • T_Ninv_T (numpy.ndarray) – Complex matrix product \(\mathbf{T}^\dagger\mathbf{N}^{-1}\mathbf{T}\). Used to construct Sigma which is used to solve for the maximum a posteriori 3D Fourier space vector and large spectral scale model amplitudes. T_Ninv_T must be a Hermitian, positive-definite matrix.

  • dbar (numpy.ndarray) – Noise weighted representation of the data vector (signal + noise) of visibilities in model (u, v, eta) space computed as \(\bar{\boldsymbol{d}}=\mathbf{T}^\dagger\mathbf{N}^{-1}\mathbf{d}\).

  • k_vals (numpy.ndarray) – Mean of each \(k\) bin.

  • k_cube_voxels_in_bin (list) – List containing sublists for each \(k\) bin. Each sublist contains the flattened 3D \(k\)-cube index of all \(\vec{k}\) that fall within a given \(k\) bin.

  • nuv (int) – Number of model uv-plane pixels per frequency channel in the EoR model.

  • neta (int) – Number of line-of-sight (frequency axis) Fourier modes.

  • nf (int) – Number of frequency channels.

  • nq (int) – Number of quadratic modes in the large spectral scale model.

  • Ninv (numpy.ndarray) – Noise covariance matrix.

  • d_Ninv_d (numpy.ndarray) – Matrix-vector product \(\boldsymbol{d}^\dagger\mathbf{N}^{-1}\boldsymbol{d}\).

  • ps_box_size_ra_Mpc (float) – Right ascension axis extent of the cosmological volume in Mpc from which the power spectrum is estimated.

  • ps_box_size_dec_Mpc (float) – Declination axis extent of the cosmological volume in Mpc from which the power spectrum is estimated.

  • ps_box_size_para_Mpc (float) – Line-of-sight extent of the cosmological volume in Mpc from which the power spectrum is estimated.

  • include_instrumental_effects (bool, optional) – Forward model an instrument. Defaults to True.

  • log_priors (bool, optional) – Assume priors on power spectrum coefficients are in log_10 units. Defaults to False.

  • uprior_inds (numpy.ndarray, optional) – Boolean 1D array that is True for any \(k\) bins using a uniform prior. False entries use a log-uniform prior. Defaults to None (all \(k\) bins use a log-uniform prior).

  • use_LWM_Gaussian_prior (bool, optional) – Use a Gaussian prior on the large spectral scale model (NOT IMPLEMENTED). Otherwise, use a uniform prior. Defaults to False.

  • inverse_LW_power (float, optional) – Inverse prior on the large spectral scale model modes. Defaults to 1e-16 (infinite variance on the large spectral scale model modes).

  • dimensionless_PS (bool, optional) – Fit for the dimensionless power spectrum, \(\Delta^2(k)\) (True), or the power spectrum, \(P(k)\) (False). Defaults to True.

  • block_T_Ninv_T (list, optional) – Block diagonal representation of T_Ninv_T. Only used if include_instrumental_effects is False. Defaults to None.

  • intrinsic_noise_fitting (bool, optional) – Fit for the amplitude of the noise in the data instead of using the noise covariance matrix, Ninv. Defaults to False.

  • fit_for_spectral_model_parameters (bool, optional) – Fit for the large spectral scale model parameter values. Defaults to False.

  • pl_max (float, optional) – Maximum brightness temperature spectral index when fitting for the optimal large spectral scale model spectral indices. Defaults to None.

  • pl_grid_spacing (float, optional) – Grid spacing for the power law spectral index axis when fitting for the large spectral scale model parameter values. Defaults to None.

  • use_shg (bool, optional) – Use the subharmonic grid. Defaults to False.

  • return_Sigma (bool, optional) – Break and return the matrix Sigma = T_Ninv_T + PhiI, where PhiI is the inverse prior covariance matrix. Defaults to False.

  • rank (int, optional) – MPI rank. Defaults to 0.

  • use_gpu (bool, optional) – Use GPUs for the Cholesky decomposition of Sigma. Otherwise, use CPUs (inadvisable due to potential inaccuracy of CPU matrix inversion). Defaults to True.

  • verbose (bool, optional) – Verbose output. Defaults to False.

  • print_rate (int, optional) – Number of iterations between print statements. Defaults to 100.

  • debug (bool, optional) – Execute break statements for debugging. Defaults to False.

  • print_debug (bool, optional) – Print debug related messages that are more detailed than verbose output. Defaults to False.

add_power_to_diagonals(T_Ninv_T_block, PhiI_block)

Add a matrix and a vector reshaped as a diagonal matrix.

Used only if T_Ninv_T has been constructed as a block-diagonal matrix, i.e. if ignoring instrumental effects.

Parameters:
  • T_Ninv_T_block (numpy.ndarray) – Single block matrix from block_T_Ninv_T.

  • PhiI_block (numpy.ndarray) – Vector of the estimated inverse variance of each \(k\) voxel present in T_Ninv_T_block.

Returns:

matrix_sum – Sum of T_Ninv_T_block and a diagonal matrix constructed from PhiI_block.

Return type:

numpy.ndarray

calc_PowerI(x)

Calculate an estimate of the inverse variance of the model.

Given a power spectrum, this function returns an estimate of the inverse variance of each voxel in the (u, v, eta) model which is used as a prior in the posterior probability calculation.

Parameters:

x (array_like) – Input power spectrum amplitudes per \(k\) bin with length nDims.

Returns:

PhiI – Vector with estimates of the inverse variance of each model (u, v, eta) voxel.

Return type:

numpy.ndarray

calc_SigmaI_dbar(Sigma, dbar, x_for_error_checking=[])

Solves the linear system Sigma * a = dbar.

Solved by calculating the Cholesky decomposition (if self.use_gpu is True) of the matrix Sigma = T_Ninv_T + PhiI. If not using GPUs to perform the matrix inversion, scipy.linalg.inv will be used. This is not desired as the CPU based scipy.linalg.inv function does not always return the “true” matrix inverse for the matrices used here.

Parameters:
  • Sigma (numpy.ndarray) – Complex matrix Sigma = T_Ninv_T + PhiI.

  • dbar (numpy.ndarray) – Noise weighted representation of the data (signal + noise) vector of visibilities in model (u, v, eta) space.

  • x_for_error_checking (array_like) – Input power spectrum amplitudes per \(k\) bin with length nDims used for error checking of the matrix inversion. Defaults to an empty list (no error checking).

Returns:

  • SigmaI_dbar (numpy.ndarray) – Complex array of maximum a posteriori (u, v, eta) and Large Spectral Scale Model (LSSM) amplitudes. Used to compute the model data vector via m = T * SigmaI_dbar.

  • logdet_Magma_Sigma (numpy.ndarray) – Natural logarith of the determinant of Sigma. Only returned if self.use_gpu is True.

calc_SigmaI_dbar_wrapper(x, T_Ninv_T, dbar, block_T_Ninv_T=[])

Wrapper for calc_SigmaI_dbar.

Constructs Sigma and calculates SigmaI_dbar and other quantities derived from SigmaI_dbar. The actual solve step to obtain SigmaI_dbar is performed in calc_SigmaI_sbar.

Parameters:
  • x (array_like) – Input power spectrum amplitudes per \(k\) bin with length nDims.

  • T_Ninv_T (numpy.ndarray) – Complex matrix product T.conjugate().T * Ninv * T.

  • dbar (numpy.ndarray) – Noise weighted representation of the data (signal + noise) vector of visibilities in model (u, v, eta) space.

  • block_T_Ninv_T (list) – Block diagonal representation of T_Ninv_T. Only used if ignoring instrumental effects. Defaults to an empty list.

Returns:

  • SigmaI_dbar (numpy.ndarray) – Complex array of maximum a posteriori (u, v, eta) and Large Spectral Scale Model (LSSM) amplitudes. Used to compute the model data vector via m = T * SigmaI_dbar.

  • dbarSigmaIdbar (numpy.ndarray) – Complex array product dbar * SigmaI_dbar.

  • PhiI (numpy.ndarray) – Vector with estimates of the inverse variance of each model (u, v, eta) voxel.

  • logSigmaDet (numpy.ndarray) – Natural logarithm of the determinant of Sigma = T_Ninv_T + PhiI.

calc_Sigma_block_diagonals(T_Ninv_T, PhiI)

Constructs Sigma using block diagonal components like block_T_Ninv_T.

Used only if T_Ninv_T has been constructed as a block-diagonal matrix, i.e. if ignoring instrumental effects.

Parameters:
  • T_Ninv_T (numpy.ndarray) – Complex matrix product T.conjugate().T * Ninv * T.

  • PhiI (numpy.ndarray) – Vector with estimates of the inverse variance of each model (u, v, eta) voxel.

Returns:

Sigma_block_diagonals – Array containing a block diagonal representation of Sigma = T_Ninv_T + PhiI. Each block is an entry along the zeroth axis of Sigma_block_diagonals. This allows Sigma to be inverted numerically block by block as opposed to all at once.

Return type:

numpy.ndarray

calc_physical_dimensionless_power_spectral_normalisation(i_bin)

Dimensionless power spectrum normalization in physical units.

This normalization will return PowerI, an estimate for one over the variance of a, in units of 1 / (mK^2 sr^2 Hz^2).

Parameters:

i_bin (int) – Spherically averaged \(k\) bin index.

Returns:

dmps_norm – Dimensionless power spectrum normalization with units of 1 / (sr**2 Hz**2).

Return type:

float

posterior_probability(x, block_T_Ninv_T=[])

Computes the posterior probability for a given power spectrum sample.

Parameters:
  • x (array_like) – Input power spectrum amplitudes per \(k\) bin with length nDims.

  • block_T_Ninv_T (list) – Block diagonal representation of T_Ninv_T. Only used if ignoring instrumental effects. Defaults to an empty list.

Returns:

  • MargLogL (float) – Posterior probability of x.

  • phi – Only used if sampling with PolyChord.

PriorC

class bayeseor.posterior.PriorC(priors_min_max)

Prior class for MultiNest and PolyChord sampler compatibility.

Parameters:

float (sequence of) – Prior [min, max] for each k bin as a a sequence, e.g. [[min1, max1], [min2, max2], …].