Power spectra (ps_eor.pspec)#

High-level Cartesian and spherical power-spectrum estimation.

Typical use: build an estimator for a cube, then ask it for the power spectrum at the dimensionality you want:

ps_gen = PowerSpectraBuilder().get(cube, z=9.1)   # -> PowerSpectraCart
ps3d = ps_gen.get_ps3d(kbins, cube)               # -> SphericalPowerSpectra
ps3d.plot()
ps3d.save_to_txt('ps3d.txt')

get_ps3d returns spherically averaged Delta^2(k); get_ps2d the cylindrical P(k_per, k_par); get_ps the per-frequency angular P(k_per); get_variance the per-frequency variance. Each returns a result object (SphericalPowerSpectra, CylindricalPowerSpectra, SpatialPowerSpectra, Variance) that knows how to plot, save, and combine with others (+ - * / propagate the error).

Units: results are stored in Kelvin^2; .get() and .plot() convert to mK^2 (or mK, the field amplitude) on request. Estimators do not subtract a noise bias – pass a noise cube to get_ps3d_with_noise for that. k is in h cMpc^-1 throughout.

class ps_eor.pspec.EorBin(name, freqs, freqs_fg, M=None)[source]#

Bases: object

Frequency range used for foreground fitting and power estimation.

get_slice(data_cube)[source]#

Return the cube slice used for power estimation.

get_slice_fg(data_cube)[source]#

Return the cube slice used for foreground fitting.

class ps_eor.pspec.EorBinList(freqs=None)[source]#

Bases: object

Named collection of frequency-bin definitions.

add_freq(name, fmhz_start, fmhz_end, fmhz_fg_start=None, fmhz_fg_end=None)[source]#

Add power and foreground-fitting frequency ranges in MHz.

Parameters:
  • name (str) – Name of the frequency bin window

  • fmhz_start (float) – Starting frequency, in MHz

  • fmhz_end (float) – Ending frequency, in MHz

  • fmhz_fg_start (float, optional) – Starting fg frequency, in MHz

  • fmhz_fg_end (float, optional) – Ending fg frequency, in MHz

get(name, freqs=None)[source]#

Resolve a named bin against an observing frequency grid.

Parameters:

name (str) – Name of the frequency bin window

Returns:

The frequency bin window

Return type:

EoRwindow

get_all(freqs=None)[source]#

Yield all resolved frequency bins.

get_all_names()[source]#

Return the bin names in insertion order.

save(filename)[source]#

Save the frequency-bin definitions as CSV.

static load(filename)[source]#

Load frequency-bin definitions from CSV.

Each non-comment row is name, fmhz_start, fmhz_end, fmhz_fg_start, fmhz_fg_end; a malformed row raises a ValueError naming the line.

Returns:

the parsed bin list.

Return type:

EorBinList

class ps_eor.pspec.MultiNightsPowerSpectraGenerator(ps_gen)[source]#

Bases: object

Apply one power-spectrum estimator to a multi-night cube.

get_variance(multi_cube, ft_nights=False, fill_gaps=True)[source]#

Estimate variance independently for each night or night mode.

get_ps2d_kpar(multi_cube, ft_nights=False)[source]#

Estimate cylindrical power and average over transverse modes.

get_ps2d_kper(multi_cube, ft_nights=False)[source]#

Estimate cylindrical power and average over line-of-sight modes.

get_ps3d(kbins, multi_cube, ft_nights=False)[source]#

Estimate spherical power independently for each night or night mode.

get_ps3d_with_noise(kbins, multi_cube, multi_cube_noise, ft_nights=False)[source]#

Estimate spherical power using paired data and noise cubes.

class ps_eor.pspec.PowerSpectraConfig(el=None, window_fct='hann', ft_method='nudft', ps2d_pos_only=True)[source]#

Bases: SimpleConfig

Estimator settings shared by the Cartesian and spherical estimators.

The knobs a user most often changes:

  • umin / umax / du – the |uv| range (in wavelengths) kept for the estimate, and the bin width for the k_per / angular binning.

  • window_fct – apodization for the frequency -> delay transform.

  • kbins_kmax / kbins_n – default spherical k-binning.

  • filter_wedge_theta / filter_kpar_min – optionally exclude the foreground wedge / low-k_par modes before spherical averaging.

  • weights_by_default – whether weighted='default' applies the uv weights (see BasePowerSpectra.get_weights()).

Load a saved config with load(); pass one to PowerSpectraBuilder, or override individual keys per call as keyword arguments to PowerSpectraBuilder.get().

class ps_eor.pspec.PowerSpectraBuilder(ps_config=None, eor_bin_list=None)[source]#

Bases: object

Factory for power-spectrum estimators.

Holds a PowerSpectraConfig (and optionally a named EorBinList); get() then returns the estimator matching a given cube and frequency bin. The usual entry point to this module.

get(cube, eor_bin_name=None, z=None, fmhz_range=None, **kargs)[source]#

Build the estimator for cube and a chosen frequency bin.

The frequency bin (which channels enter the estimate, setting the redshift and the line-of-sight comoving depth) is chosen by at most one of eor_bin_name / z / fmhz_range.

Parameters:
  • cube – a CartDataCube or SphDataCube.

  • eor_bin_name (str) – a bin name from the builder’s EorBinList.

  • z (float) – redshift; selects a +/-5 MHz window centred on it.

  • fmhz_range (tuple) – an explicit (fmin, fmax) in MHz.

  • **kargs – config keys overridden for this call.

Returns:

the estimator matching the cube type (giving no bin selector uses the full band).

Return type:

PowerSpectraCart or PowerSpectraSph

class ps_eor.pspec.BasePowerSpectra(eor_bin, ps_config, primary_beam)[source]#

Bases: object

Shared implementation for power-spectrum estimators.

get_weights(data_cube, weighted='default', delay_transform=False)[source]#

uv weights applied when estimating power from data_cube.

When on, the cube’s own weight cube is used (inverse-variance / uv coverage), and the configured wedge / k_par_min filter zeroes the excluded delay modes. Most users leave weighted='default'; it is exposed because every get_* accepts the same flag.

Parameters:
  • data_cube – the cube whose weights are read.

  • weighted'default' follows the config’s weights_by_default; True / False force weighting on / off.

  • delay_transform (bool) – return weights on the delay grid rather than the frequency grid.

Returns:

the per-mode weights, or None when unweighted.

Return type:

ndarray or None

get_all(kbins, data_cube, **kargs)[source]#

The angular, cylindrical, spherical and variance spectra of one cube.

Parameters:
  • kbins (n_k+1) – spherical k-bin edges, in h cMpc^-1.

  • data_cube – cube to estimate from.

Returns:

with .ps, .ps2d, .ps3d and .variance (see the corresponding get_* methods).

Return type:

PowerSpectraProducts

get_ps2d_norm(data_cube)[source]#

Normalization factor for 2D (spacial/frequency) PS

get_coherence_ps2d(ft_cube1, ft_cube2, cross_square=True, weighted='default')[source]#

Cylindrical coherence of two cubes: how correlated they are per (k_per, k_par) cell, in [0, 1] (dimensionless).

Formed from the cross and auto spectra as cross**2 / (P1 * P2) (cross_square=True) or cross / sqrt(P1 * P2).

Returns:

the coherence; err is zero (not propagated).

Return type:

CylindricalPowerSpectra

get_coherence_ps(ft_cube1, ft_cube2, weighted='default')[source]#

Per-frequency angular coherence of two cubes, cross**2 / (P1 * P2), in [0, 1].

Returns:

the coherence; err is zero (not propagated).

Return type:

SpatialPowerSpectra

get_coherence_variance(ft_cube1, ft_cube2, weighted='default')[source]#

Per-frequency coherence of two cubes from their variances, cross**2 / (var1 * var2), in [0, 1].

Returns:

the coherence; err is zero (not propagated).

Return type:

Variance

get_coherence_ps3d(kbins, ft_cube1, ft_cube2, cross_square=True, weighted='default')[source]#

Spherical coherence of two cubes per k-bin, cross**2 / (P1 * P2) (cross_square=True) or cross / sqrt(P1 * P2), in [0, 1].

Returns:

the coherence; err is zero (not propagated).

Return type:

SphericalPowerSpectra

get_coherence_ps3d_from_sum_diff(kbins, ft_cube1, ft_cube2, ft_cube_sum, ft_cube_diff)[source]#

Spherical coherence of two cubes using pre-summed / differenced cubes.

The cross power is taken as P(sum) - P(diff) and normalised by the autos of ft_cube1 / ft_cube2. Use when the sum and difference cubes (e.g. from interleaved data splits) are already formed.

Returns:

the coherence; err is zero (not propagated).

Return type:

SphericalPowerSpectra

get_coherence_ps2d_from_sum_diff(ft_cube1, ft_cube2, ft_cube_sum, ft_cube_diff)[source]#

Cylindrical coherence of two cubes using pre-summed / differenced cubes; cross power P(sum) - P(diff) normalised by the autos (see get_coherence_ps3d_from_sum_diff()).

Returns:

the coherence; err is zero (not propagated).

Return type:

CylindricalPowerSpectra

get_coherence_variance_from_sum_diff(ft_cube1, ft_cube2, ft_cube_sum, ft_cube_diff)[source]#

Per-frequency coherence of two cubes using pre-summed / differenced cubes; cross variance var(sum) - var(diff) normalised by the autos (see get_coherence_ps3d_from_sum_diff()).

Returns:

the coherence; err is zero (not propagated).

Return type:

Variance

get_cross_ps2d_from_sum_diff(ft_cube_sum, ft_cube_diff)[source]#

Cylindrical cross power from a summed and a differenced cube, P(sum) - P(diff) (errors added in quadrature).

Returns:

the cross power.

Return type:

CylindricalPowerSpectra

get_cross_variance_from_sum_diff(ft_cube_sum, ft_cube_diff)[source]#

Per-frequency cross variance from a summed and a differenced cube, var(sum) - var(diff) (errors added in quadrature).

Returns:

the cross variance.

Return type:

Variance

get_cross_ps3d_from_sum_diff(kbins, ft_cube_sum, ft_cube_diff)[source]#

Spherical cross power from a summed and a differenced cube, P(sum) - P(diff) (errors added in quadrature).

Parameters:
  • kbins (n_k+1) – spherical k-bin edges, in h cMpc^-1.

  • ft_cube_sum (CartDataCube) – the summed / differenced cubes.

  • ft_cube_diff (CartDataCube) – the summed / differenced cubes.

Returns:

the cross power.

Return type:

SphericalPowerSpectra

delay_transform(data_cube)[source]#

Fourier-transform the cube along frequency into delay (line-of-sight) space – the first step of every cylindrical and spherical estimate.

Applies the config window function and uv weights, and (if configured) blanks the foreground-wedge / low-k_par modes. A NoiseStdCube is handled specially (propagating its per-mode noise).

Returns:

(delay, dft_cube) – the delays (in seconds) and the delay-transformed cube.

Return type:

tuple

class ps_eor.pspec.PowerSpectraCart(eor_bin, ps_config, primary_beam)[source]#

Bases: BasePowerSpectra

Power-spectrum estimator for gridded uv (Cartesian) cubes.

Built by PowerSpectraBuilder.get() for a CartDataCube. Estimates power at three dimensionalities from the same cube – get_ps3d() (spherical), get_ps2d() (cylindrical), get_ps() (per-frequency angular) – plus get_variance(), or all of them at once with get_all(). Each get_* has a get_cross_* twin for the cross power of two cubes. All accept weighted (see get_weights()) and do not subtract a noise bias.

get_all(kbins, data_cube, weighted='default')[source]#

All four spectra in one pass, sharing the delay transform.

Cheaper than calling the individual estimators separately.

Parameters:
  • kbins (n_k+1) – spherical k-bin edges, in h cMpc^-1.

  • data_cube (CartDataCube) – cube to estimate from.

  • weighted – uv weighting, see get_weights().

Returns:

with .ps (SpatialPowerSpectra), .ps2d (CylindricalPowerSpectra), .ps3d (SphericalPowerSpectra) and .variance (Variance).

Return type:

PowerSpectraProducts

get_ps2d(data_cube, weighted='default')[source]#

Cylindrical power spectrum P(k_per, k_par).

Delay-transforms the cube and bins |V|^2 by baseline length (k_per) and delay (k_par), keeping only positive delays if ps2d_pos_only. weighted selects the uv weighting (see get_weights()).

Returns:

shape (n_kpar, n_kper), in K^2 h^-3 cMpc^3.

Return type:

CylindricalPowerSpectra

get_cross_ps2d(data_cube1, data_cube2, weighted='default', return_complex=False)[source]#

Cylindrical cross power of two cubes.

Cross power is unbiased by noise independent between the two cubes (e.g. even/odd time splits). return_complex keeps the imaginary part (otherwise the real part is taken).

Returns:

the cross power, shape (n_kpar, n_kper), in K^2 h^-3 cMpc^3; errors combine the two auto powers.

Return type:

CylindricalPowerSpectra

get_ps(data_cube, weighted='default')[source]#

Per-frequency angular power spectrum P(k_per).

The transverse power in each frequency channel (no delay transform), so its frequency structure is preserved – useful for inspecting spectral behaviour along the line of sight. weighted selects the uv weighting (see get_weights()).

Returns:

shape (n_freqs, n_el), in K^2 h^-2 cMpc^2.

Return type:

SpatialPowerSpectra

get_cl(data_cube, weighted='default')[source]#

Angular spectrum C_l – get_ps() divided by X^2 (the angular-to-comoving factor).

Returns:

C_l per frequency (cl=True), in K^2 sr.

Return type:

SpatialPowerSpectra

get_cross_ps(data_cube1, data_cube2, weighted='default', return_complex=False)[source]#

Per-frequency angular cross power of two cubes (see get_cross_ps2d() on cross power; return_complex keeps the imaginary part).

Returns:

shape (n_freqs, n_el), in K^2 h^-2 cMpc^2.

Return type:

SpatialPowerSpectra

get_variance(data_cube, weighted='default')[source]#

Variance of the sky signal per frequency channel.

The k_per-integral of the angular power (get_cl()), i.e. the total transverse variance in each channel.

Returns:

shape (n_freqs,), in K^2 (plot in mK^2).

Return type:

Variance

get_cross_variance(data_cube1, data_cube2, weighted='default')[source]#

Per-frequency cross variance of two cubes (see get_cross_ps2d()).

Returns:

shape (n_freqs,), in K^2.

Return type:

Variance

get_ps3d(kbins, data_cube, weighted='default')[source]#

Spherically averaged dimensionless power Delta^2(k).

Delay-transforms the cube, restricts to the config uv range [umin, umax], and averages the power in the spherical shells given by kbins. This is the usual final data product.

Parameters:
  • kbins (n_k+1) – spherical k-bin edges, in h cMpc^-1.

  • data_cube (CartDataCube) – cube to estimate from.

  • weighted – uv weighting, see get_weights().

Returns:

Delta^2(k), shape (n_k,), in K^2 (read out in mK^2 via SphericalPowerSpectra.get()). Not noise-debiased – see get_ps3d_with_noise().

Return type:

SphericalPowerSpectra

get_ps3d_with_noise(kbins, ft_cube, noise_cube, weighted='default')[source]#

Noise-debiased spherical power: Delta^2(data) - Delta^2(noise).

Parameters:
  • kbins (n_k+1) – spherical k-bin edges, in h cMpc^-1.

  • ft_cube (CartDataCube) – the data cube.

  • noise_cube (CartDataCube) – a noise realization or expectation, same geometry as ft_cube.

  • weighted – uv weighting, see get_weights().

Returns:

the debiased Delta^2(k); the two estimates’ errors add in quadrature.

Return type:

SphericalPowerSpectra

get_cross_ps3d(kbins, data_cube1, data_cube2, weighted='default', return_complex=False)[source]#

Spherically averaged cross power of two cubes (see get_cross_ps2d() on cross power; return_complex keeps the imaginary part).

Returns:

Delta^2(k), shape (n_k,), in K^2.

Return type:

SphericalPowerSpectra

class ps_eor.pspec.PowerSpectraSph(eor_bin, ps_config, primary_beam)[source]#

Bases: BasePowerSpectra

Power-spectrum estimator for spherical-harmonic cubes.

get_ps2d(alm_cube)[source]#

Cylindrical power P(k_per, k_par) from a spherical-harmonic cube.

The spherical-cube counterpart of PowerSpectraCart.get_ps2d(); errors are the analytic sample variance sqrt(2 / ((2l+1) f_sky)).

Returns:

shape (n_kpar, n_kper), in K^2 h^-3 cMpc^3.

Return type:

CylindricalPowerSpectra

get_ps(alm_cube)[source]#

Per-frequency angular power P(k_per) from a spherical-harmonic cube (counterpart of PowerSpectraCart.get_ps(); analytic sample-variance errors).

Returns:

shape (n_freqs, n_el), in K^2 h^-2 cMpc^2.

Return type:

SpatialPowerSpectra

get_ps3d(kbins, alm_cube)[source]#

Spherically averaged Delta^2(k) from a spherical-harmonic cube (counterpart of PowerSpectraCart.get_ps3d()).

Parameters:
  • kbins (n_k+1) – spherical k-bin edges, in h cMpc^-1.

  • alm_cube (SphDataCube) – the spherical-harmonic cube.

Returns:

Delta^2(k), shape (n_k,), in K^2.

Return type:

SphericalPowerSpectra

get_ps3d_with_noise(kbins, ft_cube, noise_cube)[source]#

Noise-debiased spherical power, Delta^2(data) - Delta^2(noise) (counterpart of PowerSpectraCart.get_ps3d_with_noise()).

Returns:

the debiased Delta^2(k).

Return type:

SphericalPowerSpectra

get_cl(alm_cube)[source]#

Angular spectrum C_l from a spherical-harmonic cube (counterpart of PowerSpectraCart.get_cl()).

Returns:

C_l per frequency (cl=True), in K^2 sr.

Return type:

SpatialPowerSpectra

get_variance(alm_cube)[source]#

Per-frequency variance from a spherical-harmonic cube (counterpart of PowerSpectraCart.get_variance()).

Returns:

shape (n_freqs,), in K^2.

Return type:

Variance

class ps_eor.pspec.PowerSpectraMath[source]#

Bases: object

Element-wise arithmetic for the result containers, with error propagation.

Lets you combine spectra directly – a - b to subtract a model or a noise estimate, 0.5 * a to rescale, a / b for a ratio – returning the same type with the 1-sigma error propagated (added in quadrature for +/-, the standard ratio formula for /). Operands must share the same binning.

class ps_eor.pspec.Variance(var, var_err, freqs, var_w=1)[source]#

Bases: PowerSpectraMath

Sky-signal variance per frequency channel.

data, err

variance and its 1-sigma error, in K^2.

Type:

n_freqs

freqs#

channel frequencies, in Hz.

Type:

n_freqs

plot() (converts to mK^2), freq_binning() to rebin in frequency; supports arithmetic via PowerSpectraMath.

freq_binning(df)[source]#

Average the result into frequency bins of width df.

plot(ax=None, df=None, nsigma=1, title=None, mkelvin=True, **kargs)[source]#

Plot the variance vs frequency on a log-y axis.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • df (float) – if set, average into frequency bins of this width first.

  • nsigma (int) – error-bar width.

  • mkelvin (bool) – plot in mK^2 (else K^2).

  • **kargs – forwarded to the matplotlib errorbar call.

class ps_eor.pspec.SpatialPowerSpectra(ps, ps_err, freqs, el, k_per, cl=False, n_eff=None, ps_w=1)[source]#

Bases: PowerSpectraMath

Angular power P(k_per) per frequency channel (from PowerSpectraCart.get_ps()).

data, err

P(k_per) and its 1-sigma error, per channel, in K^2 h^-2 cMpc^2.

Type:

n_freqs, n_el

freqs#

channel frequencies (Hz). el (n_el): angular multipoles.

Type:

n_freqs

k_per#

transverse wavenumber (h cMpc^-1). n_eff: effective mode count per cell. cl: True if holding C_l instead of P(k_per).

Type:

n_el

plot() (frequency vs k_per image), plot_kper(), save_to_txt() / load(); arithmetic via PowerSpectraMath.

plot(ax=None, title=None, k_only=True, fill_gap=True, text=None, log_norm=True, l_lambda=False, normalize=False, imaginary_part=False, **kargs)[source]#

Plot P(k_per) as a frequency vs k_per image.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • k_only (bool) – label the transverse axis in k_per only (no twin l axis).

  • fill_gap (bool) – interpolate across missing frequency channels.

  • log_norm (bool) – colour scale in log.

  • **kargs – forwarded to the matplotlib imshow call.

plot_kper(ax=None, nsigma=0, fill_std=False, normalize=False, mkelvin=True, kelvin_square=True, weighted=True, l_lambda=False, **kargs)[source]#

Plot P(k_per), averaged over frequency, as a 1-D curve vs k_per.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • nsigma (int) – error-bar / band width (0 for the line only).

  • weighted (bool) – weight the frequency average by the uv weights.

  • mkelvin – units, see SphericalPowerSpectra.get().

  • kelvin_square – units, see SphericalPowerSpectra.get().

  • **kargs – forwarded to the matplotlib call.

save_to_txt(filename)[source]#

Write a plain-text table (one row per freq x k_per cell: freq, k_per, baseline, P and error in K^2 h^-2 cMpc^2, N_eff). See load().

static load(filename, z=None)[source]#

Read a save_to_txt() table back into a SpatialPowerSpectra; z recovers the multipoles from k_per when the file lists only k_per.

class ps_eor.pspec.CylindricalPowerSpectra(ps2d, ps2d_err, delay, el, k_per, k_par, n_eff=None, ps2d_w=1)[source]#

Bases: PowerSpectraMath

Cylindrical power P(k_per, k_par) (from PowerSpectraCart.get_ps2d()).

The 2-D power in transverse (k_per) and line-of-sight (k_par) scale – the plane where the foreground “wedge” is diagnosed.

data, err

P and its 1-sigma error, in K^2 h^-3 cMpc^3.

Type:

n_kpar, n_kper

k_per#

wavenumbers (h cMpc^-1). delay (n_kpar): the line-of-sight delay (us). el (n_kper): angular multipoles.

Type:

n_kper), k_par (n_kpar

n_eff#

effective mode count per cell.

plot() (2-D image), plot_kpar() / plot_kper() (1-D cuts), reduce_region() (summarise a k-space box), save_to_txt() / load(); arithmetic via PowerSpectraMath.

reduce_region(kpar=None, kper=None, reducer='median')[source]#

Reduce power within a rectangular region of cylindrical k-space.

Bounds are (lower, upper) with an inclusive lower edge and an exclusive upper edge. Either edge may be None. Non-finite values are ignored and an empty region returns nan.

plot(ax=None, title=None, k_only=True, log_norm=True, colorbar=True, log_axis=False, ax_cb=None, text=None, dimensionless=False, wedge_lines=None, z=None, imaginary_part=False, **kargs)[source]#

Plot P(k_per, k_par) as a 2-D image – the standard wedge diagnostic.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • k_only (bool) – label the k axes only (no twin delay / baseline axes).

  • log_norm (bool) – colour scale in log.

  • colorbar (bool) – add a colour bar.

  • log_axis (bool) – log-scale the k axes.

  • wedge_lines (list) – angles (deg) at which to overlay wedge lines.

  • dimensionless (bool) – show Delta^2 instead of P.

  • **kargs – forwarded to imshow.

plot_kpar(ax=None, nsigma=0, fill_std=False, delay=False, weighted=True, **kargs)[source]#

Plot P(k_par), averaged over k_per, as a 1-D curve vs k_par.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • nsigma (int) – error-bar / band width (0 for the line only).

  • delay (bool) – use delay (us) on the x-axis instead of k_par.

  • weighted (bool) – weight the k_per average by the uv weights.

  • **kargs – forwarded to the matplotlib call.

plot_kper(ax=None, nsigma=0, fill_std=False, normalize=False, weighted=True, **kargs)[source]#

Plot P(k_per), averaged over k_par, as a 1-D curve vs k_per.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • nsigma (int) – error-bar / band width (0 for the line only).

  • normalize (bool) – show the dimensionless Delta^2(k_per) instead of P.

  • weighted (bool) – weight the k_par average by the uv weights.

  • **kargs – forwarded to the matplotlib call.

save_to_txt(filename)[source]#

Write a plain-text table (one row per k_par x k_per cell: k_par, k_per, delay, baseline, P and error in K^2 h^-3 cMpc^3, N_eff). See load().

static load(filename, z=None)[source]#

Read a save_to_txt() table back into a CylindricalPowerSpectra; z recovers delay/multipoles from k_par/k_per when only those are listed.

class ps_eor.pspec.SphericalPowerSpectra(ps3d, ps3d_err, k_bins, k_mean, ps3d_q16=None, ps3d_q84=None, n_eff=None, k_std=None)[source]#

Bases: PowerSpectraMath

Spherically averaged dimensionless power Delta^2(k) (from PowerSpectraCart.get_ps3d()).

The usual final 21-cm data product: one number per k-bin.

data, err

Delta^2(k) and its 1-sigma error, in K^2. When 16/84 posterior quantiles are provided, err is the half 68% width.

Type:

n_k

k_mean, k_std

mean and spread of |k| within each bin (h cMpc^-1).

Type:

n_k

k_bins#

bin edges. n_eff (n_k): effective mode count per bin.

Type:

n_k+1

get() returns the spectrum in the requested units (mK^2, or mK for the amplitude); get_upper() an n-sigma upper limit; plot(); save_to_txt() / load_from_txt(). Subtract a noise/model estimate with - (see PowerSpectraMath).

get(mkelvin=True, kelvin_square=False)[source]#

The spectrum and its error in the requested units.

Parameters:
  • mkelvin (bool) – scale K -> mK.

  • kelvin_square (bool) – return the power Delta^2 (True) or its square root, the amplitude Delta (False, the default).

Returns:

(data, err) arrays. The default returns Delta(k) in mK.

Return type:

tuple

get_upper(nsigma=2, mkelvin=True, kelvin_square=False)[source]#

The data + nsigma*err upper limit (units per get()).

Returns:

the upper limit per k-bin.

Return type:

ndarray

plot(ax=None, nsigma=2, marker='+', mkelvin=True, kelvin_square=True, title=None, fill_std=False, kerr_as_kbins=False, imaginary_part=False, **kargs)[source]#

Plot Delta^2(k) vs k on log-log axes.

Parameters:
  • ax – matplotlib Axes to draw on (a new figure if None).

  • nsigma (int) – error-bar / band width; 0 draws the line only.

  • mkelvin – units, see get().

  • kelvin_square – units, see get().

  • fill_std (bool) – shade a +/-nsigma band instead of drawing error bars.

  • **kargs – forwarded to the matplotlib plot / errorbar call.

save_to_txt(filename)[source]#

Write a plain-text table (one row per k-bin: k edges, k_mean, k_std, Delta^2 and error in mK^2, N_eff). Round-trips through load_from_txt().

static load_from_txt(filename)[source]#

Read a save_to_txt() table back into a SphericalPowerSpectra (mK^2 in the file -> K^2 internally). Also accepts a 3-column k_mean/Delta^2/error file.

class ps_eor.pspec.FourPanelPsResults(ps_gen, kbins, figsize=(10, 8))[source]#

Bases: object

Build a 2x2 comparison figure: variance, spherical Delta^2(k), P(k_par) and P(k_per).

Overlay one or more cubes / stackers with the add_* methods (each takes a legend label), then call done() and savefig().

add_cube(cube, label, ps_gen=None, **kargs)[source]#

Estimate and overlay one cube’s four spectra, labelled label.

add_cube_ps_diff(cube1, cube2, label, ps_gen=None, **kargs)[source]#

Overlay the difference of two cubes (cube1 - cube2; noise-debiased for the Delta^2(k) panel).

add_ps_stacker(ps_stacker, label, **kargs)[source]#

Overlay a PsStacker’s four spectra, with credible bands.

done(ncol_legend=3)[source]#

Add the shared legend and tighten the layout. Call before savefig().

savefig(filename, **kargs)[source]#

Save the figure to filename (run done() first).

class ps_eor.pspec.ThreePanelPsResults(ps_gen, kbins, figsize=(10, 5), norm_factor=1, n_cols=1, dpi=100)[source]#

Bases: object

Build a 3-row comparison figure (Delta^2(k), P(k_par), P(k_per)), with an optional grid of columns and a shared norm_factor.

Overlay cubes / stackers with the add_* methods, then done() and savefig(). col selects the column to draw into.

add_cube(cube, label, ps_gen=None, col=0, **kargs)[source]#

Estimate and overlay one cube’s three spectra in column col.

add_cube_ps_diff(cube1, cube2, label, ps_gen=None, col=0, **kargs)[source]#

Overlay the difference of two cubes (noise-debiased Delta^2(k)) in column col.

add_ps_stacker(ps_stacker, label, col=0, **kargs)[source]#

Overlay a PsStacker’s three spectra in column col.

done(ncol_legend=3, col_legend=0, **kargs)[source]#

Add the shared legend and tighten the layout. Call before savefig().

savefig(filename, **kargs)[source]#

Save the figure to filename (run done() first).

class ps_eor.pspec.MultiNight2DPowerSpectra(ps, err, x, xlabel, y, ylabel)[source]#

Bases: PowerSpectraMath

A quantity binned on two generic axes (e.g. a spectrum vs observing night), carrying its own x / y values and axis labels.

data, err

the 2-D values and 1-sigma error. x, y: the axis coordinates (y may be string category labels). xlabel, ylabel: display labels.

Plot with plot(); supports arithmetic via PowerSpectraMath.

plot(ax=None, log_norm=True, colorbar=True, **kargs)[source]#

Draw the 2-D map as an image.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • log_norm (bool) – colour scale in log.

  • colorbar (bool) – add a colour bar.

  • **kargs – forwarded to imshow.

class ps_eor.pspec.PowerSpectraMC(all_ps)[source]#

Bases: object

A stack of power-spectrum realizations, summarised by median and credible bands.

Built from a list of same-binning spectra (typically one per posterior sample). Exposes the sample median (med, used as data by the dimensioned subclasses), robust std, and the 16/84 and 2.5/97.5 percentiles (q16 / q84 / q2_5 / q97_5) for 68% / 95% bands. Arithmetic combines two stacks by resampling, propagating the spread.

Use the dimensioned subclasses – SphericalPowerSpectraMC, CylindricalPowerSpectraMC, SpatialPowerSpectraMC, VarianceMC – which add the matching plot / save / load.

get_data(mkelvin=False)[source]#

The median and credible-band percentiles in the requested units.

Returns:

(median, q16, q84, q2_5, q97_5, unit) – the sample median and the 16/84/2.5/97.5 percentiles, in K (or mK if mkelvin).

Return type:

tuple

class ps_eor.pspec.SphericalPowerSpectraMC(ps3d_all)[source]#

Bases: PowerSpectraMC, SphericalPowerSpectra

Posterior stack of spherical Delta^2(k): a SphericalPowerSpectra whose data is the sample median, carrying 68% / 95% credible bands (see PowerSpectraMC). Returned by PsStacker.get_ps3d().

get_upper(nsigma=2, mkelvin=True, kelvin_square=False)[source]#

The posterior upper limit from the credible-band percentiles: the 84th (nsigma=1) or 97.5th (nsigma=2) percentile.

Returns:

the upper limit per k-bin, units per get().

Return type:

ndarray

plot(ax=None, show68=True, show95=True, marker='+', mkelvin=True, title=None, kerr_as_kbins=False, **kargs)[source]#

Plot the median Delta^2(k) with shaded 68% / 95% credible bands.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • show68 (bool) – shade the 68% / 95% band.

  • show95 (bool) – shade the 68% / 95% band.

  • mkelvin (bool) – plot in mK^2 (else K^2).

  • **kargs – forwarded to the matplotlib plot call.

save(filename)[source]#

Save every realization to an HDF5 file (and a .txt summary alongside). Round-trips through load().

save_to_txt(filename)[source]#

Write a text table: median Delta^2 and the 2.5/16/84/97.5 percentiles (all in mK^2) per k-bin.

static load(filename)[source]#

Load a save() HDF5 file back into a SphericalPowerSpectraMC.

class ps_eor.pspec.CylindricalPowerSpectraMC(all_ps2d)[source]#

Bases: PowerSpectraMC, CylindricalPowerSpectra

Posterior stack of cylindrical P(k_per, k_par): a CylindricalPowerSpectra at the sample median, with 68% / 95% bands (see PowerSpectraMC). Returned by PsStacker.get_ps2d().

plot_kpar(ax=None, show68=True, show95=True, delay=False, weighted=True, mkelvin=False, **kargs)[source]#

Plot the median P(k_par) (averaged over k_per) with 68% / 95% bands.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • show68 (bool) – shade the 68% / 95% band.

  • show95 (bool) – shade the 68% / 95% band.

  • delay (bool) – use delay (us) on the x-axis instead of k_par.

  • weighted (bool) – weight the k_per average by the uv weights.

  • **kargs – forwarded to the matplotlib call.

plot_kper(ax=None, show68=True, show95=True, normalize=False, weighted=True, mkelvin=False, **kargs)[source]#

Plot the median P(k_per) (averaged over k_par) with 68% / 95% bands.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • show68 (bool) – shade the 68% / 95% band.

  • show95 (bool) – shade the 68% / 95% band.

  • normalize (bool) – show the dimensionless Delta^2(k_per) instead of P.

  • weighted (bool) – weight the k_par average by the uv weights.

  • **kargs – forwarded to the matplotlib call.

save(filename)[source]#

Save every realization to an HDF5 file. Round-trips through load().

static load(filename)[source]#

Load a save() HDF5 file back into a CylindricalPowerSpectraMC.

class ps_eor.pspec.VarianceMC(all_var)[source]#

Bases: PowerSpectraMC, Variance

Posterior stack of per-frequency variance: a Variance at the sample median, with 68% / 95% bands (see PowerSpectraMC). Returned by PsStacker.get_variance().

plot(ax=None, show68=True, show95=True, mkelvin=True, title=None, **kargs)[source]#

Plot the median variance vs frequency with 68% / 95% credible bands.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • show68 (bool) – shade the 68% / 95% band.

  • show95 (bool) – shade the 68% / 95% band.

  • mkelvin (bool) – plot in mK^2 (else K^2).

  • **kargs – forwarded to the matplotlib call.

save(filename)[source]#

Save every realization to an HDF5 file. Round-trips through load().

static load(filename)[source]#

Load a save() HDF5 file back into a VarianceMC.

class ps_eor.pspec.SpatialPowerSpectraMC(all_var)[source]#

Bases: PowerSpectraMC, SpatialPowerSpectra

Posterior stack of per-frequency angular P(k_per): a SpatialPowerSpectra at the sample median, with 68% / 95% bands (see PowerSpectraMC). Returned by PsStacker.get_ps().

plot_kper(ax=None, show68=True, show95=True, normalize=False, weighted=True, l_lambda=False, mkelvin=True, **kargs)[source]#

Plot the median P(k_per) (or C_l) averaged over frequency, with 68% / 95% credible bands.

Parameters:
  • ax – matplotlib Axes (new figure if None).

  • show68 (bool) – shade the 68% / 95% band.

  • show95 (bool) – shade the 68% / 95% band.

  • normalize (bool) – show the dimensionless form (Delta^2 or l(l+1)C_l).

  • weighted (bool) – weight the frequency average by the uv weights.

  • **kargs – forwarded to the matplotlib call.

save(filename)[source]#

Save every realization to an HDF5 file. Round-trips through load().

static load(filename)[source]#

Load a save() HDF5 file back into a SpatialPowerSpectraMC.

class ps_eor.pspec.PowerSpectraProducts(ps: SpatialPowerSpectra, ps2d: CylindricalPowerSpectra, ps3d: SphericalPowerSpectra, variance: Variance)[source]#

Bases: object

Power-spectrum products computed from one cube.

class ps_eor.pspec.PsStacker(ps_gen, kbins)[source]#

Bases: object

Accumulate power spectra from repeated cube realizations.

add(cube)[source]#

Estimate and retain all supported spectra for one cube.

get_ps()[source]#

Per-frequency angular power P(k_perp), as a SpatialPowerSpectraMC over the stacked realizations (axes: frequency MHz x k_perp h cMpc^-1; units K^2 h^-2 cMpc^2). No noise-bias subtraction: it is the power of whatever cubes were added.

get_ps2d()[source]#

Cylindrical power P(k_perp, k_par), as a CylindricalPowerSpectraMC over the stacked realizations (axes k_perp, k_par in h cMpc^-1; units K^2 h^-3 cMpc^3). No noise-bias subtraction.

get_variance()[source]#

Per-frequency variance, as a VarianceMC over the stacked realizations (in K^2). The spread across the stack is the posterior uncertainty.

get_ps3d()[source]#

Spherically averaged dimensionless power Delta^2(k), as a SphericalPowerSpectraMC over the stacked realizations (binned in the PsStacker’s kbins, k in h cMpc^-1; units mK^2). No noise-bias subtraction. The spread across the stack is the posterior uncertainty.

save(dir_path, name)[source]#

Save all accumulated result families.

static load(dir_path, name)[source]#

Load a previously saved stack.