UVPSpec Class

UVPSpec objects contain power spectra and associated metadata.

class hera_pspec.UVPSpec[source]

An object for storing power spectra generated by hera_pspec and a file-format for its data and meta-data.

antnums_to_bl(antnums)[source]

Convert tuple of antenna numbers to baseline integer.

Parameters

antnums (tuple) – Tuple containing integer antenna numbers for a baseline. Ex. (ant1, ant2)

Returns

bl – Baseline integer.

Return type

i6 integer

antnums_to_blpair(antnums)[source]

Convert nested tuple of antenna numbers to baseline-pair integer.

Parameters

antnums (tuple) – nested tuple containing integer antenna numbers for a baseline-pair. Ex. ((ant1, ant2), (ant3, ant4))

Returns

blpair – baseline-pair integer

Return type

i12 integer

average_spectra(blpair_groups=None, time_avg=False, blpair_weights=None, error_field=None, error_weights=None, inplace=True, add_to_history='')[source]

Average power spectra across the baseline-pair-time axis, weighted by each spectrum’s integration time.

This is an “incoherent” average, in the sense that this averages power spectra, rather than visibility data. The nsample_array and integration_array will be updated to reflect the averaging.

In the case of averaging across baseline pairs, the resultant averaged spectrum is assigned to the zeroth blpair in the group. In the case of time averaging, the time and LST arrays are assigned to the mean of the averaging window.

Note that this is designed to be separate from spherical binning in k: here we are not connecting k_perp modes to k_para modes. However, if blpairs holds groups of iso baseline separation, then this is equivalent to cylindrical binning in 3D k-space.

If you want help constructing baseline-pair groups from baseline groups, see self.get_blpair_groups_from_bl_groups.

Parameters

blpair_groups (list) – List of list of baseline-pair group tuples or integers. All power spectra in a baseline-pair group are averaged together. If a baseline-pair exists in more than one group, a warning is raised.

Examples:

blpair_groups = [ [((1, 2), (1, 2)), ((2, 3), (2, 3))],
                  [((4, 6), (4, 6))]]

blpair_groups = [ [1002001002, 2003002003], [4006004006] ]
time_avgbool, optional

If True, average power spectra across the time axis. Default: False.

blpair_weightslist, optional

List of float or int weights dictating the relative weight of each baseline-pair when performing the average.

This is useful for bootstrapping. This should have the same shape as blpair_groups if specified. The weights are automatically normalized within each baseline-pair group. Default: None (all baseline pairs have unity weights).

error_field: string or list, optional

If errorbars have been entered into stats_array, will do a weighted sum to shrink the error bars down to the size of the averaged data_array. Error_field strings be keys of stats_array. If list, does this for every specified key. Every stats_array key that is not specified is thrown out of the new averaged object.

error_weights: string, optional

error_weights specify which kind of errors we use for weights during averaging power spectra. The weights are defined as $w_i = 1/ sigma_i^2$, where $sigma_i$ is taken from the relevant field of stats_array.

If error_weight is set to None, which means we just use the integration time as weights. If error_weights is specified, then it also gets appended to error_field as a list. Default: None.

inplacebool, optional

If True, edit data in self, else make a copy and return. Default: True.

add_to_historystr, optional

Added text to add to file history.

Notes

Currently, every baseline-pair in a blpair group must have the same Ntpairs, unless time_avg=True. Future versions may support baseline-pair averaging of heterogeneous time arrays. This includes the scenario of repeated blpairs (e.g. in bootstrapping), which will return multiple copies of their time_array.

bl_to_antnums(bl)[source]

Convert baseline (antenna-pair) integer to nested tuple of antenna numbers.

Parameters

bl (i6 int) – baseline integer ID

Returns

antnums – tuple containing baseline antenna numbers. Ex. (ant1, ant2)

Return type

tuple

blpair_to_antnums(blpair)[source]

Convert baseline-pair integer to nested tuple of antenna numbers.

Parameters

blpair (i12 int) – baseline-pair integer ID

Returns

antnums – Nested tuple containing baseline-pair antenna numbers. Ex. ((ant1, ant2), (ant3, ant4))

Return type

tuple

blpair_to_indices(blpair)[source]

Convert a baseline-pair nested tuple ((ant1, ant2), (ant3, ant4)) or a baseline-pair integer into indices to index the blpairts axis of data_array.

Parameters

blpair (tuples or ints) – Nested tuple or blpair i12 integer, Ex. ((1, 2), (3, 4)), or list of blpairs.

check(just_meta=False)[source]

Run checks to make sure metadata and data arrays are properly defined and in the right format. Raises AssertionErrors if checks fail.

Parameters

just_meta (bool, optional) – If True, only check metadata. Default: False.

compute_scalar(spw, polpair, num_steps=1000, little_h=True, noise_scalar=False)[source]

Compute power spectrum normalization scalar given an adopted cosmology and a beam model. See pspecbeam.PSpecBeamBase.compute_pspec_scalar for details.

Parameters
  • spw (integer) – Spectral window selection.

  • polpair (tuple or int or str) – Which polarization pair to calculate the scalar for.

  • num_steps (int, optional) – Number of integration bins along frequency in computing scalar. Default: 1000.

  • little_h (bool, optional) – Whether to use h^-1 units or not. Default: True.

  • noise_scalar (bool, optional) – If True calculate noise pspec scalar, else calculate normal pspec scalar. See pspecbeam.py for difference between normal scalar and noise scalar. Default: False.

Returns

scalar – Power spectrum normalization scalar.

Return type

float

convert_to_deltasq(inplace=True)[source]

Convert from P(k) to Delta^2(k) by multiplying by k^3 / (2pi^2).

The units of the output is therefore the current units (self.units) times (h^3) Mpc^-3, where the h^3 is only included if h^-3 is in self.norm_units

Parameters
  • inplace (bool, optional) – If True edit and overwrite arrays in self, else make a copy of self and return. Default: True.

  • Notes

  • ——

  • Alters data_array, stats_array, and cov_array if present,

  • as well as metadata like norm_units

fold_spectra()[source]

Average bandpowers from matching positive and negative delay bins onto a purely positive delay axis. Negative delay bins are still populated, but are filled with zero power. This is an in-place operation.

Will only work if self.folded == False, i.e. data is currently unfolded across negative and positive delay. Because this averages the data, the nsample array is multiplied by a factor of 2.

WARNING: This operation cannot be undone.

generate_noise_spectra(spw, polpair, Tsys, blpairs=None, little_h=True, form='Pk', num_steps=2000, component='real', scalar=None)[source]

Generate the expected RMS noise power spectrum given a selection of spectral window, system temp. [K], and polarization. This estimate is constructed as:

P_N = scalar * (Tsys * 1e3)^2 / (integration_time) / sqrt(Nincoherent)

[mK^2 h^-3 Mpc^3]

where scalar is the cosmological and beam scalar (i.e. X2Y * Omega_eff) calculated from pspecbeam with noise_scalar = True, integration_time is in seconds and comes from self.integration_array and Nincoherent is the number of incoherent averaging samples and comes from self.nsample_array. If component is real or imag, P_N is divided by an additional factor of sqrt(2).

If the polarizations specified are pseudo Stokes pol (I, Q, U or V) then an extra factor of 2 is divided.

If form is DelSq then a factor of |k|^3 / (2pi^2) is multiplied.

If real is True, a factor of sqrt(2) is divided to account for discarding imaginary noise component.

For more details, see hera_pspec.noise.Sensitivity.calc_P_N, and Cheng et al. 2018.

The default units of P_N are mK^2 h^-3 Mpc^3 (if little_h=True and form=’Pk’), but is different if those kwargs are altered.

Parameters
  • spw (int) – Spectral window index to generate noise curve for.

  • polpair (tuple or int or str) – Polarization-pair selection in form of tuple (e.g. (‘I’, ‘Q’)) or polpair int. Strings are expanded into polarization pairs, e.g. ‘XX’ becomes (‘XX,’XX’).

  • Tsys (dictionary, float or array) – System temperature in Kelvin for each blpair. Key is blpair-integer, value is Tsys float or ndarray. If fed as an ndarray, shape=(Ntpairs,)

  • blpairs (list) – List of unique blair tuples or i12 integers to calculate noise spectrum for default is to calculate for baseline pairs.

  • little_h (bool, optional) – Whether to have cosmological length units be h^-1 Mpc or Mpc Default: h^-1 Mpc

  • form (str, optional) – Form of power spectra, whetehr P(k) or Delta^2(k). Valid options are ‘Pk’ or ‘DelSq’. Default: ‘Pk’.

  • num_steps (int, optional) – Number of frequency bins to use in integrating power spectrum scalar in pspecbeam. Default: 2000.

  • component (str) – Options=[‘real’, ‘imag’, ‘abs’]. If component is real or imag, divide by an extra factor of sqrt(2).

  • scalar (float) – Optional noise scalar used to convert from Tsys to cosmological units output from pspecbeam.compute_pspec_scalar(…,noise_scalar=True) Default is None. If None provided, will calculate scalar in function.

Returns

P_N – Dictionary containing blpair integers as keys and float ndarrays of noise power spectra as values, with ndarrays having shape (Ntimes, Ndlys).

Return type

dict

get_ENU_bl_vecs()[source]

Return baseline vector array in TOPO (ENU) frame in meters, with matched ordering of self.bl_vecs.

Returns

blvecs – Array of baseline vectors.

Return type

array_like

get_all_keys()[source]

Returns a list of all possible tuple keys in the data_array, in the format:

(spectral window, baseline-pair, polarization-string)

get_blpair_blvecs(use_second_bl=False)[source]

For each baseline-pair, get the baseline vector in ENU frame of the first baseline in the pair. Ordering matches self.get_blpairs()

Parameters

use_second_bl (bool) – If True, return the vector of the second bl in the blpair.

Returns

vector in ENU frame of the baseline

Return type

ndarray shape=(Nblpairs,)

get_blpair_groups_from_bl_groups(blgroups, only_pairs_in_bls=False)[source]

Get baseline pair matches from a list of baseline groups.

Parameters
  • blgroups (list) – List of baseline groups, which themselves are lists of baseline tuples or baseline i6 integers. Ex: [ [(1, 2), (2, 3), (3, 4)], [(1, 4), (2, 5)] ]

  • only_pairs_in_bls (bool, optional) – If True, select only baseline-pairs whose first _and_ second baseline are both found in each baseline group. Default: False.

Returns

blpair_groups – List of blpair groups, which themselves are lists of blpair integers.

Return type

list

get_blpair_seps()[source]

For each baseline-pair, get the average baseline separation in ENU frame in meters between its bl1 and bl2. Ordering matches self.get_blpairs()

Returns

blp_avg_sep – shape=(Nbltpairs,)

Return type

float ndarray

get_blpairs()[source]

Returns a list of all blpair tuples in the data_array.

get_cov(key, component='real', omit_flags=False)[source]

Slice into covariance array with a specified data key in the format (spw, ((ant1, ant2),(ant3, ant4)), (pol1, pol2))

or

(spw, blpair-integer, pol-integer)

where spw is the spectral window integer, ant1 etc. are integers, and pol is either a tuple of polarization strings (ex. (‘XX’, ‘YY’)) or integers (ex. (-5,-5)).

Parameters
  • key (tuple) – Contains the baseline-pair, spw, polpair keys

  • component (str) – “real” or “imag”. Indicating which cov_array the function calls.

  • omit_flags (bool, optional) – If True, remove time integrations (or spectra) that came from visibility data that were completely flagged across the spectral window (i.e. integration == 0).

Returns

data – Shape (Ntimes, Ndlys, Ndlys)

Return type

complex ndarray

get_data(key, omit_flags=False)[source]

Slice into data_array with a specified data key in the format

(spw, ((ant1, ant2), (ant3, ant4)), (pol1, pol2))

or

(spw, blpair-integer, polpair)

where spw is the spectral window integer, ant1 etc. are integers, and polpair is either a tuple of polarization strings/ints or a polarizatio-pair integer.

Parameters
  • key (tuple) – Contains the baseline-pair, spw, polpair keys

  • omit_flags (bool, optional) – If True, remove time integrations (or spectra) that came from visibility data that were completely flagged across the spectral window (i.e. integration == 0).

Returns

data – Shape (Ntimes, Ndlys)

Return type

complex ndarray

get_dlys(key)[source]

Get array of delays given a spectral window selection.

Parameters

key (int, or tuple with integer) – Spectral window selection

Returns

dlys

Return type

float ndarray, contains delays of pspectra in nanosec, given spw

get_exact_window_functions(ftbeam=None, spw_array=None, verbose=False, inplace=True, add_to_history='', x_orientation=None)[source]

Obtain the exact window functions corresponding to UVPSpec object.

Depending on the format of the UVPSpec object given as input, the shape of the output array will be different: (Nblpts, Ndlys, Nkperp, Nkpara, Npols). Weights can be applied to the different baseline pairs to facilitate spherical average later.

Parameters
  • ftbeam (str or FTBeam object, optional) – Definition of the beam Fourier transform to be used. Options include;

    • Root name of the file to use, without the polarisation

    Ex : FT_beam_HERA_dipole (+ path) - ‘’ for computation from beam simulations (slow) - FTBeam object. Make sure the polarisation and bandwidths are compatible with uvp.

  • spw_array (list of ints, optional) – Spectral window indices. If None, the window functions are computed on all the uvp.spw_ranges, successively. Default: None.

  • verbose (bool, optional) – If True, print progress, warnings and debugging info to stdout.

  • inplace (bool, optional) – If True (default value), the UVPspec attribute window_function_array is filled with the values computed in the function, and window_function_kperp_bins and window_function_kpara_bins array are added as attributes. If False, returns kperp_bins, kpara_bins and window functions computed. Automatically set to False if blpair_groups is not None (that is, if the window functions are not computed on all the blpairs).

  • add_to_history (str, optional) – Added text to add to file history.

  • x_orientation (str, optional) – Orientation in cardinal direction east or north of X dipole. Default keeps polarization in X and Y basis. Used to convert polstr to polnum and conversely. Default: None.

get_integrations(key, omit_flags=False)[source]

Slice into integration_array with a specified data key in the format:

(spw, ((ant1, ant2), (ant3, ant4)), (pol1, pol2))

or

(spw, blpair-integer, polpair-integer)

where spw is the spectral window integer, ant1 etc. are integers, and polpair is either a polarization pair tuple or integer.

Parameters
  • key (tuple) – Contains the baseline-pair, spw, polpair keys

  • omit_flags (bool, optional) – If True, remove time integrations (or spectra) that came from visibility data that were completely flagged across the spectral window (i.e. integration == 0).

Returns

data – Has shape (Ntimes,)

Return type

float ndarray

get_kparas(spw, little_h=True)[source]

Get radial (parallel) cosmological wavevectors for power spectra given an adopted cosmology and a spw selection.

Parameters
  • spw (int) – Choice of spectral window

  • little_h (boolean, optional) – Whether to have cosmological length units be h^-1 Mpc or Mpc Default: h^-1 Mpc

  • Returns (k_perp, k_para)

  • ——-

  • k_para (float ndarray) – Radial wave-vectors, shape=(Ndlys given spw)

get_kperps(spw, little_h=True)[source]

Get transverse (perpendicular) cosmological wavevector for each baseline-pair given an adopted cosmology and a spw selection.

Parameters
  • spw (int, choice of spectral window)

  • little_h (boolean, optional) – Whether to have cosmological length units be h^-1 Mpc or Mpc Default: h^-1 Mpc

Returns

k_perp

Return type

float ndarray, transverse wave-vectors, shape=(Nblpairs,)

get_nsamples(key, omit_flags=False)[source]

Slice into nsample_array with a specified data key in the format

(spw, ((ant1, ant2), (ant3, ant4)), (pol1, pol2))

or

(spw, blpair-integer, polpair-integer)

where spw is the spectral window integer, ant1 etc. are integers, and polpair is either a polarization pair tuple or integer.

Parameters
  • key (tuple) – Contains the baseline-pair, spw, polpair keys

  • omit_flags (bool, optional) – If True, remove time integrations (or spectra) that came from visibility data that were completely flagged across the spectral window (i.e. integration == 0).

Returns

data

Return type

float ndarray with shape (Ntimes,)

get_polpairs()[source]

Returns a list of all unique polarization pairs in the data_array.

get_r_params()[source]

Return an r_params dictionary.

In a PSpecData object, the r_params are stored as a dictionary with one entry per baseline.

In a UVPSpec object, the dictionary is compressed so that a single r_param entry corresponds to multiple baselines and is stored as a JSON format string.

This function reads the compressed string and returns the dictionary with the correct following fields and structure.

Returns

r_params – Dictionary of r_params for this object.

Return type

dict

get_spw_ranges(spws=None)[source]

Get the frequency range and Nfreqs of each spectral window in the object.

Parameters

spws (list of spw integers, optional) – Default is all spws present in data.

Returns

spw_ranges – Tuples: (freq_start, freq_end, Nfreqs, Ndlys) in Hz. Contains start, stop and Nfreqs and Ndlys [Hz] of each spectral window. To turn this into the frequency array of the spectral window use spw_freqs = np.linspace(*spw_range[:3], endpoint=False)

Return type

list of len-4 tuples

get_stats(stat, key, omit_flags=False)[source]

Returns a statistic from the stats_array dictionary.

Parameters
  • stat (str) – The statistic to return.

  • key (tuple) – A spw-blpair-polpair key parseable by UVPSpec.key_to_indices.

  • omit_flags (bool, optional) – If True, remove time integrations (or spectra) that came from visibility data that were completely flagged across the spectral window (i.e. integration == 0).

Returns

statistic – A 2D (Ntimes, Ndlys) complex array holding desired bandpower statistic.

Return type

array_like

get_wgts(key, omit_flags=False)[source]

Slice into wgt_array with a specified data key in the format

(spw, ((ant1, ant2), (ant3, ant4)), (pol1, pol2))

or

(spw, blpair-integer, polpair-integer)

where spw is the spectral window integer, ant1 etc. are integers, and polpair is either a polarization pair tuple or integer.

Parameters
  • key (tuple) – Contains the baseline-pair, spw, polpair keys

  • omit_flags (bool, optional) – If True, remove time integrations (or spectra) that came from visibility data that were completely flagged across the spectral window (i.e. integration == 0).

Returns

wgts – Has shape (Ntimes, Nfreqs, 2), where the last axis holds [wgt_1, wgt_2] in that order.

Return type

float ndarray

get_window_function(key, omit_flags=False)[source]

Slice into window_function array with a specified data key in the format (spw, ((ant1, ant2),(ant3, ant4)), (pol1, pol2))

or

(spw, blpair-integer, pol-integer)

where spw is the spectral window integer, ant1 etc. are integers, and pol is either a tuple of polarization strings (ex. (‘XX’, ‘YY’)) or integers (ex. (-5,-5)).

Parameters
  • key (tuple) – Contains the baseline-pair, spw, polpair keys

  • omit_flags (bool, optional) – If True, remove time integrations (or spectra) that came from visibility data that were completely flagged across the spectral window (i.e. integration == 0).

Returns

data – Shape (Ntimes, Ndlys, Ndlys)

Return type

complex ndarray

key_to_indices(key, omit_flags=False)[source]

Convert a data key into relevant slice arrays. A data key takes the form

(spw_integer, ((ant1, ant2), (ant3, ant4)), (pol1, pol2))

or

(spw, blpair-integer, pol-integer)

where spw is the spectral window integer, ant1 etc. are integers, and polpair is either a polarization-pair integer or tuple.

The key can also be a dictionary in the form:

key = {
  'spw' : spw_integer,
  'blpair' : ((ant1, ant2), (ant3, ant4))
  'polpair' : (pol1, pol2)
  }

and it will parse the dictionary for you.

Parameters
  • key (tuple) – Baseline-pair, spw, and pol-pair key.

  • omit_flags (bool, optional) – If True, remove time integrations (or spectra) that came from visibility data that were completely flagged across the spectral window (i.e. integration == 0).

Returns

  • spw (int) – Spectral window index.

  • blpairts (list) – List of integers to apply along blpairts axis.

  • polpair (int) – Polarization pair index.

polpair_to_indices(polpair)[source]

Map a polarization-pair integer or tuple to its index in polpair_array.

Parameters

polpair ((list of) int or tuple or str) – Polarization-pair integer or tuple of polarization strings/ints.

Alternatively, if a single string is given, will assume that the specified polarization is to be cross-correlated withitself, e.g. ‘XX’ implies (‘XX’, ‘XX’).

Returns

indices – Index of polpair in polpair_array.

Return type

int

read_from_group(grp, just_meta=False, spws=None, bls=None, blpairs=None, times=None, lsts=None, polpairs=None, only_pairs_in_bls=False)[source]

Clear current UVPSpec object and load in data from specified HDF5 group.

Parameters
  • grp (HDF5 group) – HDF5 group to load data from.

  • just_meta (bool, optional) – If True, read-in metadata but ignore data, wgts and integration arrays. Default: False.

  • spws (list of tuple, optional) – List of spectral window integers to select. Default: None (loads all channels).

  • bls (list of i6 baseline integers or baseline tuples) – Select all baseline-pairs whose first _or_ second baseline are in the list. This changes if only_pairs_in_bls == True. Example tuple: (2, 3). Default: None (loads all bl pairs).

  • blpairs (list of baseline-pair tuples or integers) – List of baseline pairs to keep. If bls is also fed, this list is concatenated onto the baseline-pair list constructed from from the bls selection.

  • times (float ndarray) – Times from the time_avg_array to keep.

  • lsts (float ndarray) – lsts from the lst_avg_array to keep.

  • polpairs (list of tuple or int or str) – List of polarization-pair integers, tuples, or strings to keep. Strings are expanded into polarization pairs, and are not treated as individual polarizations, e.g. ‘XX’ becomes (‘XX,’XX’).

  • only_pairs_in_bls (bool, optional) – If True, keep only baseline-pairs whose first _and_ second baseline are found in the ‘bls’ list. Default: False.

read_hdf5(filepath, just_meta=False, spws=None, bls=None, blpairs=None, times=None, lsts=None, polpairs=None, only_pairs_in_bls=False)[source]

Clear current UVPSpec object and load in data from an HDF5 file.

Parameters
  • filepath (str) – Path to HDF5 file to read from.

  • just_meta (bool, optional) – If True, read-in metadata but ignore data, wgts and integration arrays. Default: False.

  • spws (list of tuple, optional) – List of spectral window integers to select. Default: None (loads all channels).

  • bls (list of i6 baseline integers or baseline tuples) – Select all baseline-pairs whose first _or_ second baseline are in the list. This changes if only_pairs_in_bls == True. Example tuple: (2, 3). Default: None (loads all bl pairs).

  • blpairs (list of baseline-pair tuples or integers) – List of baseline pairs to keep. If bls is also fed, this list is concatenated onto the baseline-pair list constructed from from the bls selection.

  • times (float ndarray) – Times from the time_avg_array to keep.

  • lsts (float ndarray) – lsts from the lst_avg_array to keep.

  • polpairs (list of polpair ints or tuples or str) – List of polarization-pair integers or tuples to keep.

    Strings are expanded into polarization pairs, and are not treated as individual polarizations, e.g. ‘XX’ becomes (‘XX,’XX’).

  • only_pairs_in_bls (bool, optional) – If True, keep only baseline-pairs whose first _and_ second baseline are found in the ‘bls’ list. Default: False.

select(spws=None, bls=None, only_pairs_in_bls=True, blpairs=None, times=None, lsts=None, polpairs=None, inplace=True, run_check=True)[source]

Select function for selecting out certain slices of the data.

Parameters
  • spws (list, optional) – List of spectral window integers to select. If None, all will be selected. Default: None.

  • bls (list of i6 baseline integers or baseline tuples, optional) – Example: (2,3). Select all baseline-pairs whose first _or_ second baseline are in bls list. This changes if only_pairs_in_bls == True. If None, all baselines are selected (subject to other options). Default: None.

  • only_pairs_in_bls (bool, optional) – If True, keep only baseline-pairs whose first _and_ second baseline are found in bls list. Default: True.

  • blpairs (list of baseline-pair tuples or integers, optional) – List of baseline-pairs to keep. If bls is also fed, this list is concatenated onto the baseline-pair list constructed from the bls selection. If None, all valid baseline pairs are selected. Default: None.

  • times (array_like, optional) – Float ndarray of times from the time_avg_array to select. If None, all times are kept. Default: None.

  • lsts (array_like, optional) – Float ndarray of lsts from the lst_avg_array to select. If None, all lsts are kept. Default: None.

  • polpairs (list of tuple or int or str, optional) – List of polarization-pairs to keep. If None, all polarizations are kept. Default: None.

    (Strings are expanded into polarization pairs, and are not treated as individual polarizations, e.g. ‘XX’ becomes (‘XX,’XX’).)

  • inplace (bool, optional) – If True, edit and overwrite arrays in self, else make a copy of self and return. Default: True.

  • run_check (bool, optional) – If True, run uvp.check() after selection

Returns

uvp – If inplace=False, return a new UVPSpec object containing only the selected data.

Return type

UVPSpec, optional

set_cosmology(new_cosmo, overwrite=False, new_beam=None, verbose=True)[source]

Set the cosmology for this UVPSpec object by passing an instance of conversions.Cosmo_Conversions and assigning it to self.cosmo, in addition to re-computing power spectrum normalization quantities like self.scalar_array. Because this will attempt to recompute the scalar_array, the beam-related metadata (OmegaP, OmegaPP and beam_freqs) must exist in self.

If they do not, or if you’d like to overwrite them with a new beam, you can pass a UVBeam object or path to a beam file via the new_beam kwarg.

If self.cosmo already exists then you are attempting to overwrite the currently adopted cosmology, which will only proceed if overwrite is True.

If overwrite == True, then this module will overwrite self.cosmo. It will also recompute the power spectrum normalization scalar (using pspecbeam) and will overwrite the values in self.scalar_array. It will then propagate these re-normalization changes into the data_array by multiplying the data by (new_scalar_array / old_scalar_array).

Parameters
  • new_cosmo (Cosmo_Conversions or dict) – Cosmo_Conversions instance or cosmological parameter dictionary. The new cosmology you want to adopt for this UVPSpec object.

  • overwrite (bool, optional) – If True, overwrite self.cosmo if it already exists. Default: False.

  • new_beam (PSpecBeamBase sublcass or str) – pspecbeam.PSpecBeamUV object or path to beam file. The new beam you want to adopt for this UVPSpec object.

  • verbose (bool, optional) – If True, report feedback to stdout. Default: True.

set_stats(stat, key, statistic)[source]

Set a statistic waterfall (Ntimes, Ndlys) in the stats_array given the selection of spw, blpairts and polpair in key.

Parameters
  • stat (string) – Name of the statistic.

  • key (tuple) – A spw-blpair-polpair key parseable by UVPSpec.key_to_indices.

  • statistic (ndarray) – 2D ndarray with statistics to set. Must conform to shape of the slice of data_array produced by the specified key.

set_stats_slice(stat, m, b, above=False, val=1e+20)[source]

For each baseline, set all delay bins in stats_array that fall above or below y = bl_len * m + b [nanosec] equal to val. Useful for downweighting foregrounds in spherical average.

Parameters
  • stat (str, name of stat in stat_array to set)

  • m (float, coefficient of bl_len [nanosec / meter])

  • b (float, offset [nanosec])

  • above (bool, if True, set stats above line, else set below line)

  • val (float, value to replace in stats_array)

spw_indices(spw)[source]

Convert a spectral window integer into a list of indices to index into the spw_array.

Parameters

spw (int or tuple) – Spectral window index, or spw tuple from get_spw_ranges(), or list of either.

Returns

spw_indices – Index array for specified spw(s) in spw_array.

Return type

array_like

spw_to_dly_indices(spw)[source]

Convert a spectral window integer into a list of indices to index into the spw_dly_array or dly_array.

Parameters

spw (int or tuple) – Spectral window index, or spw tuple from get_spw_ranges(), or list of either.

Returns

dly_indices – Index array for specified spw(s) in spw_dly_array.

Return type

array_like

spw_to_freq_indices(spw)[source]

Convert a spectral window integer into a list of indices to index into the spw_freq_array or freq_array.

Parameters

spw (int or tuple) – Spectral window index, or spw tuple from get_spw_ranges(), or list of either.

Returns

freq_indices – Index array for specified spw(s) in spw_freq_array.

Return type

array_like

time_to_indices(time, blpairs=None)[source]

Convert a time [Julian Date] from self.time_avg_array to an array that indexes the elements at which it appears in that array. Can optionally take a blpair selection to further select the indices at which both that time and blpair are satisfied.

Parameters
  • time (float) – Julian Date time from self.time_avg_array, Ex: 2458042.12242

  • blpairs (tuple or int, optional) – List of blpair tuples or integers that further selects the elements at which both time and blpairs are satisfied.

Returns

indices – Contains indices at which selection is satisfied along the blpairts axis.

Return type

integer ndarray

property units

Return power spectrum units. See self.vis_units and self.norm_units.

write_hdf5(filepath, overwrite=False, run_check=True)[source]

Write a UVPSpec object to HDF5 file.

Parameters
  • filepath (str) – Filepath for output file.

  • overwrite (bool, optional) – Whether to overwrite output file if it exists. Default: False.

  • run_check (bool, optional) – Run UVPSpec validity check before writing to file. Default: True.

write_to_group(group, run_check=True)[source]

Write UVPSpec data into an HDF5 group.

Parameters
  • group (HDF5 group) – The handle of the HDF5 group that the UVPSpec object should be written into.

  • run_check (bool, optional) – Whether to run a validity check on the UVPSpec object before writing it to the HDF5 group. Default: True.