UVBeam Class

UVBeam is the main user class for beam models. It provides import and export functionality to supported file formats and can be interacted with directly.

class pyuvdata.UVBeam[source]

A class for defining a radio telescope antenna beam.

UVParameter objects

For full list see UVBeam Parameters (http://pyuvdata.readthedocs.io/en/latest/uvbeam_parameters.html). Some are always required, some are required for certain beam_types, antenna_types and pixel_coordinate_systems and others are always optional.

check(check_extra=True, run_check_acceptability=True)[source]

Check that all required parameters are set reasonably.

Check that required parameters exist and have appropriate shapes. Optionally check if the values are acceptable.

Parameters
  • check_extra – Option to check optional parameters as well as required ones. Default is True.

  • run_check_acceptability – Option to check if values in required parameters are acceptable. Default is True.

set_cs_params()[source]

Set various forms and required parameters depending on pixel_coordinate_system.

set_efield()[source]

Set beam_type to ‘efield’ and adjust required parameters.

set_power()[source]

Set beam_type to ‘power’ and adjust required parameters.

set_simple()[source]

Set antenna_type to ‘simple’ and adjust required parameters.

set_phased_array()[source]

Set antenna_type to ‘phased_array’ and adjust required parameters.

peak_normalize()[source]

Convert to peak normalization.

efield_to_pstokes(run_check=True, check_extra=True, run_check_acceptability=True, inplace=True)[source]

Convert E-field to pseudo-stokes power as done in https://arxiv.org/pdf/1802.04151.pdf.

M_ij = Tr(sigma_i J sigma_j J^*)

where sigma_i and sigma_j are Pauli matrices.

Parameters
  • run_check – Option to check for the existence and proper shapes of the required parameters after converting to power. Default is True.

  • run_check_acceptability – Option to check acceptable range of the values of required parameters after combining objects. Default is True.

  • check_extra – Option to check optional parameters as well as required ones. Default is True.

  • inplace – Option to perform the select directly on self (True, default) or return a new UVBeam object, which is a subselection of self (False).

efield_to_power(calc_cross_pols=True, keep_basis_vector=False, run_check=True, check_extra=True, run_check_acceptability=True, inplace=True)[source]

Convert E-field beam to power beam.

Parameters
  • calc_cross_pols – If True, calculate the crossed polarization beams (e.g. ‘xy’ and ‘yx’), otherwise only calculate the same polarization beams (e.g. ‘xx’ and ‘yy’). Default is True.

  • keep_basis_vector – If True, keep the directionality information and just multiply the efields for each basis vector separately (caution: this is not what is standardly meant by the power beam). Default is False.

  • run_check – Option to check for the existence and proper shapes of required parameters after converting to power. Default is True.

  • check_extra – Option to check optional parameters as well as required ones. Default is True.

  • run_check_acceptability – Option to check acceptable range of the values of required parameters after combining objects. Default is True.

  • inplace – Option to perform the select directly on self (True, default) or return a new UVBeam object, which is a subselection of self (False)

interp(az_array=None, za_array=None, freq_array=None, freq_interp_tol=1.0, polarizations=None, freq_interp_kind='linear', reuse_spline=False)[source]

Interpolate beam to given az, za locations (in radians).

Parameters
  • az_array – az values to interpolate to (same length as za_array)

  • za_array – za values to interpolate to (same length as az_array)

  • freq_array – frequency values to interpolate to

  • freq_interp_tol – float, frequency distance tolerance [Hz] of nearest neighbors. If all elements in freq_array have nearest neighbor distances within the specified tolerance then return the beam at each nearest neighbor, otherwise interpolate the beam.

  • polarizations – list of str, polarizations to interpolate if beam_type is ‘power’. Default is all polarizations in self.polarization_array.

  • freq_interp_kind – str, interpolation method across frequency. See scipy.interpolate.interp1d for details.

  • reuse_spline – Save the interpolation functions for reuse. Only applies for az_za_simple interpolation.

Returns

an array of interpolated values, shape

(Naxes_vec, Nspws, Nfeeds or Npols,

Nfreqs or freq_array.size if freq_array is passed, Npixels/(Naxis1, Naxis2) or az_array.size if az/za_arrays are passed)

an array of interpolated basis vectors (or self.basis_vector_array

if az/za_arrays are not passed), shape: (Naxes_vec, Ncomponents_vec, Npixels/(Naxis1, Naxis2) or az_array.size if az/za_arrays are passed)

to_healpix(nside=None, run_check=True, check_extra=True, run_check_acceptability=True, inplace=True)[source]

Convert beam in to healpix coordinates. The interpolation is done using the interpolation method specified in self.interpolation_function.

Note that this interpolation isn’t perfect. Interpolating an Efield beam and then converting to power gives a different result than converting to power and then interpolating at about a 5% level.

Parameters
  • nside – The nside to use for the Healpix map. If not specified, use the nside that gives the closest resolution that is higher than the input resolution.

  • run_check – Option to check for the existence and proper shapes of required parameters after converting to healpix. Default is True.

  • check_extra – Option to check optional parameters as well as required ones. Default is True.

  • run_check_acceptability – Option to check acceptable range of the values of required parameters after combining objects. Default is True.

  • inplace – Option to perform the select directly on self (True, default) or return a new UVBeam object, which is a subselection of self (False)

get_beam_area(pol='pI')[source]

Computes the integral of the beam, which has units of steradians

Pseudo-Stokes ‘pI’ (I), ‘pQ’(Q), ‘pU’(U), ‘pV’(V) beam and linear dipole ‘XX’, ‘XY’, ‘YX’ and ‘YY’ are supported. See Equations 4 and 5 of Moore et al. (2017) ApJ 836, 154 or arxiv:1502.05072 and Kohn et al. (2018) or https://arxiv.org/pdf/1802.04151.pdf for details.

Parameters

pol – polarization string, Ex. a pseudo-stokes pol ‘pI’, or a linear pol ‘XX’

Returns

omega – float, integral of the beam across the sky [steradians]

get_beam_sq_area(pol='pI')[source]

Computes the integral of the beam^2, which has units of steradians

Pseudo-Stokes ‘pI’ (I), ‘pQ’(Q), ‘pU’(U), ‘pV’(V) beam and linear dipole ‘XX’, ‘XY’, ‘YX’ and ‘YY’ are supported. See Equations 4 and 5 of Moore et al. (2017) ApJ 836, 154 or arxiv:1502.05072 for details.

Parameters

pol – polarization string, Ex. a pseudo-stokes pol ‘pI’, or a linear pol ‘XX’

Returns

omega – float, integral of the beam^2 across the sky [steradians]

select(axis1_inds=None, axis2_inds=None, pixels=None, frequencies=None, freq_chans=None, feeds=None, polarizations=None, run_check=True, check_extra=True, run_check_acceptability=True, inplace=True)[source]

Select specific image axis indices or pixels (if healpix), frequencies and feeds or polarizations (if power) to keep in the object while discarding others.

The history attribute on the object will be updated to identify the operations performed.

Parameters
  • axis1_inds – The indices along the first image axis to keep in the object. Cannot be set if pixel_coordinate_system is “healpix”.

  • axis2_inds – The indices along the second image axis to keep in the object. Cannot be set if pixel_coordinate_system is “healpix”.

  • pixels – The healpix pixels to keep in the object. Cannot be set if pixel_coordinate_system is not “healpix”.

  • frequencies – The frequencies to keep in the object.

  • freq_chans – The frequency channel numbers to keep in the object.

  • feeds – The feeds to keep in the object. Cannot be set if the beam_type is “power”.

  • polarizations – The polarizations to keep in the object. Cannot be set if the beam_type is “efield”.

  • run_check – Option to check for the existence and proper shapes of required parameters after downselecting data on this object. Default is True.

  • check_extra – Option to check optional parameters as well as required ones. Default is True.

  • run_check_acceptability – Option to check acceptable range of the values of required parameters after downselecting data on this object. Default is True.

  • inplace – Option to perform the select directly on self (True, default) or return a new UVBeam object, which is a subselection of self (False)

read_beamfits(filename, run_check=True, check_extra=True, run_check_acceptability=True)[source]

Read in data from a beamfits file.

Parameters
  • filename – The beamfits file or list of files to read from.

  • run_check – Option to check for the existence and proper shapes of required parameters after reading in the file. Default is True.

  • check_extra – Option to check optional parameters as well as required ones. Default is True.

  • run_check_acceptability – Option to check acceptable range of the values of required parameters after reading in the file. Default is True.

write_beamfits(filename, run_check=True, check_extra=True, run_check_acceptability=True, clobber=False)[source]

Write the data to a beamfits file.

Parameters
  • filename – The beamfits file to write to.

  • run_check – Option to check for the existence and proper shapes of required parameters before writing the file. Default is True.

  • check_extra – Option to check optional parameters as well as required ones. Default is True.

  • run_check_acceptability – Option to check acceptable range of the values of required parameters before writing the file. Default is True.

  • clobber – Option to overwrite the filename if the file already exists. Default is False.

read_cst_beam(filename, beam_type='power', feed_pol=None, rotate_pol=None, frequency=None, telescope_name=None, feed_name=None, feed_version=None, model_name=None, model_version=None, history=None, x_orientation=None, reference_impedance=None, extra_keywords=None, frequency_select=None, run_check=True, check_extra=True, run_check_acceptability=True)[source]

Read in data from a cst file.

Parameters
  • filename (str) – Either a settings yaml file or a cst text file or list of cst text files to read from. If a list is passed, the files are combined along the appropriate axes. Settings yaml files must include the following keywords:

    - telescope_name (str)
    - feed_name (str)
    - feed_version (str)
    - model_name (str)
    - model_version (str)
    - history (str)
    - frequencies (list(float))
    - cst text filenames (list(str)) – path relative to yaml file location
    - feed_pol (str) or (list(str))
    and they may include the following optional keywords:
    - x_orientation (str): Optional but strongly encouraged!
    - ref_imp (float): beam model reference impedance
    - sim_beam_type (str): e.g. ‘E-farfield’
    - all other fields will go into the extra_keywords attribute

    More details and an example are available in the docs (cst_settings_yaml.rst). Specifying any of the associated keywords to this function will override the values in the settings file.

  • beam_type (str) – what beam_type to read in (‘power’ or ‘efield’). Defaults to ‘power’.

  • feed_pol (str) – the feed or polarization or list of feeds or polarizations the files correspond to. Defaults to ‘x’ (meaning x for efield or xx for power beams).

  • rotate_pol (bool) – If True, assume the structure in the simulation is symmetric under 90 degree rotations about the z-axis (so that the y polarization can be constructed by rotating the x polarization or vice versa). Default: True if feed_pol is a single value or a list with all the same values in it, False if it is a list with varying values.

  • frequency (list(float)) – the frequency or list of frequencies corresponding to the filename(s). This is assumed to be in the same order as the files. If not passed, the code attempts to parse it from the filenames.

  • telescope_name (str) – the name of the telescope corresponding to the filename(s).

  • feed_name (str) – the name of the feed corresponding to the filename(s).

  • feed_version (str) – the version of the feed corresponding to the filename(s).

  • model_name (str) – the name of the model corresponding to the filename(s).

  • model_version (str) – the version of the model corresponding to the filename(s).

  • history (str) – A string detailing the history of the filename(s).

  • x_orientation (str) – Orientation of the physical dipole corresponding to what is labelled as the x polarization. Options are “east” (indicating east/west orientation) and “north” (indicating north/south orientation)

  • reference_impedance (float) – The reference impedance of the model(s).

  • extra_keywords (dict) – a dictionary containing any extra_keywords.

  • frequency_select (list(float)) – Only used if the file is a yaml file. Indicates which frequencies to include (only read in files for those frequencies)

  • run_check – Option to check for the existence and proper shapes of required parameters after reading in the file. Default is True.

  • check_extra – Option to check optional parameters as well as required ones. Default is True.

  • run_check_acceptability – Option to check acceptable range of the values of required parameters after reading in the file. Default is True.