UVCal

UVCal is the main user class for calibration solutions for interferometric data sets. It provides import and export functionality to and from the supported file formats (calfits, FHD) as well as methods for transforming the data (converting types, selecting, sorting) and can be interacted with directly.

Attributes

The attributes on UVCal hold all of the metadata and data required to work with calibration solutions for interferometric data sets. Under the hood, the attributes are implemented as properties based on pyuvdata.parameter.UVParameter objects but this is fairly transparent to users.

UVCal objects can be initialized as an empty object (as cal = UVCal()). When an empty UVCal object is initialized, it has all of these attributes defined but set to None. The attributes can be set by reading in a data file using the pyuvdata.UVCal.read_calfits() or pyuvdata.UVCal.read_fhd_cal() methods or by setting them directly on the object. Some of these attributes are required to be set to have a fully defined calibration data set while others are optional. The pyuvdata.UVCal.check() method can be called on the object to verify that all of the required attributes have been set in a consistent way.

Note that objects can be in a “metadata only” state where all of the metadata is defined but the data-like attributes (gain_array, delay_array, flag_array, quality_array) are not. The pyuvdata.UVCal.check() method will still pass for metadata only objects.

Note location type attributes (which are given in topocentric xyz coordinates) have convenience properties named the same thing with _lat_lon_alt and _lat_lon_alt_degrees appended through which you can get or set the values using latitude, longitude and altitude values in radians or degrees and meters.

Required

These parameters are required to have a sensible UVCal object and are required for most kinds of uv cal files.

Nants_data

Number of antennas that have data associated with them (i.e. length of ant_array), which may be smaller than the numberof antennas in the telescope (i.e. length of antenna_numbers).

Nants_telescope

Number of antennas in the antenna_numbers array. May be larger than the number of antennas with gains associated with them.

Nfreqs

Number of frequency channels

Njones

Number of Jones calibration parameters (Number of Jones matrix elements calculated in calibration).

Nspws

Number of spectral windows (ie non-contiguous spectral chunks).

Ntimes

Number of times with different calibrations calculated (if a calibration is calculated over a range of integrations, this gives the number of separate calibrations along the time axis).

ant_array

Array of integer antenna numbers that appear in self.gain_array, with shape (Nants_data,). This array is ordered to match the inherent ordering of the zeroth axis of self.gain_array.

antenna_names

Array of antenna names with shape (Nants_telescope,). Ordering of elements matches ordering of antenna_numbers.

antenna_numbers

Array of all integer-valued antenna numbers in the telescope with shape (Nants_telescope,). Ordering of elements matches that of antenna_names. This array is not necessarily identical to ant_array, in that this array holds all antenna numbers associated with the telescope, not just antennas with data, and has an in principle non-specific ordering.

antenna_positions

Array giving coordinates of antennas relative to telescope_location (ITRF frame), shape (Nants_telescope, 3), units meters. See the tutorial page in the documentation for an example of how to convert this to topocentric frame.

cal_style

Style of calibration. Values are sky or redundant.

cal_type

cal type parameter. Values are delay or gain.

channel_width

Width of frequency channels (Hz). If flex_spw = False and future_array_shapes=False, then it is a single value of type = float, otherwise it is an array of shape (Nfreqs,), type = float.Should not be set if future_array_shapes=True and wide_band=True.

flex_spw

Option to construct a ‘flexible spectral window’, which storesall spectral channels across the frequency axis of data_array. Allows for spectral windows of variable sizes, and channels of varying widths.

freq_array

Array of frequencies, center of the channel, shape (1, Nfreqs) or (Nfreqs,) if future_array_shapes=True, units Hz.Not required if future_array_shapes=True and wide_band=True.Should not be set if future_array_shapes=True and wide_band=True.

future_array_shapes

Flag indicating that this object is using the future array shapes.

gain_convention

The convention for applying the calibration solutions to data.Values are “divide” or “multiply”, indicating that to calibrate one should divide or multiply uncalibrated data by gains. Mathematically this indicates the alpha exponent in the equation: calibrated data = gain^alpha * uncalibrated data. A value of “divide” represents alpha=-1 and “multiply” represents alpha=1.

history

String of history, units English

integration_time

Integration time of a time bin, units seconds. If future_array_shapes=False, then it is a single value of type = float, otherwise it is an array of shape (Ntimes), type = float.

jones_array

Array of antenna polarization integers, shape (Njones). linear pols -5:-8 (jxx, jyy, jxy, jyx).circular pols -1:-4 (jrr, jll. jrl, jlr).

lst_array

Array of lsts, center of integration, shape (Ntimes), units radians

spw_array

Array of spectral window numbers, shape (Nspws).

telescope_location

Telescope location: xyz in ITRF (earth-centered frame). Can also be accessed using telescope_location_lat_lon_alt or telescope_location_lat_lon_alt_degrees properties

telescope_name

Name of telescope. e.g. HERA. String.

time_array

Array of calibration solution times, center of integration, shape (Ntimes), units Julian Date

wide_band

Option to support ‘wide-band’ calibration solutions with gains or delays that apply over a range of frequencies rather than having distinct values at each frequency. Delay type cal solutions are always ‘wide-band’ if future_array_shapes is True. If it is True several other parameters are affected: future_array_shapes is also True; the data-like arrays have a spw axis that is Nspws long rather than a frequency axis that is Nfreqs long; the freq_range parameter is required and the freq_array and channel_width parameters are not required.

x_orientation

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)

Optional

These parameters are defined by one or more file standard but are not always required. Some of them are required depending on the cal_type or cal_style (as noted below).

Nsources

Number of sources used.

baseline_range

Range of baselines used for calibration.

delay_array

Required if cal_type = “delay”. Array of delays with units of seconds. Shape: (Nants_data, 1, 1, Ntimes, Njones) or (Nants_data, Nspws, Ntimes, Njones) if future_array_shapes=True, type=float.

diffuse_model

Name of diffuse model.

extra_keywords

Any user supplied extra keywords, type=dict. Keys should be 8 character or less strings if writing to calfits files. Use the special key ‘comment’ for long multi-line string comments.

filename

List of strings containing the unique basenames (not the full path) of input files.

flag_array

Array of flags to be applied to calibrated data (logical OR of input and flag generated by calibration). True is flagged. Shape: (Nants_data, 1, Nfreqs, Ntimes, Njones) or (Nants_data, Nfreqs, Ntimes, Njones) if future_array_shapes=True and wide_band=False or (Nants_data, Nspws, Ntimes, Njones) if wide_band=True, type = bool.

flex_spw_id_array

Required if flex_spw = True and will always be required for non-wide-band objects starting in version 3.0. Maps individual channels along the frequency axis to individual spectral windows, as listed in the spw_array. Shape (Nfreqs), type = int.

freq_range

Required if cal_type=’delay’ or wide_band=True. Frequency range that solutions are valid for. If future_array_shapes is False it is a length 2 array with [start_frequency, end_frequency], otherwise it is an array of shape (Nspws, 2). Units are Hz.Should not be set if cal_type=’gain’ and wide_band=False.

gain_array

Required if cal_type = “gain”. Array of gains, shape: (Nants_data, 1, Nfreqs, Ntimes, Njones) or (Nants_data, Nfreqs, Ntimes, Njones) if future_array_shapes=True, or (Nants_data, Nspws, Ntimes, Njones) if wide_band=True, type = complex float.

gain_scale

The gain scale of the calibration, which indicates the units of the calibrated visibilities. For example, Jy or K str.

git_hash_cal

Commit hash of calibration software (from git_origin_cal) used to generate solutions.

git_origin_cal

Origin (on github for e.g) of calibration software. Url and branch.

input_flag_array

Deprecated, support will be removed in version 2.5. Array of input flags, True is flagged. shape: (Nants_data, 1, Nfreqs, Ntimes, Njones) or (Nants_data, Nfreqs, Ntimes, Njones) if future_array_shapes=True, type = bool.

observer

Name of observer who calculated solutions in this file.

quality_array

Array of qualities of calibration solutions. The shape depends on cal_type, if the cal_type is ‘gain’, the shape is: (Nants_data, 1, Nfreqs, Ntimes, Njones) or (Nants_data, Nfreqs, Ntimes, Njones) if future_array_shapes=True and wide_band=False or (Nants_data, Nspws, Ntimes, Njones) if wide_band=True, if the cal_type is ‘delay’, the shape is (Nants_data, 1, 1, Ntimes, Njones) or (Nants_data, Nspws, Ntimes, Njones) if future_array_shapes=True. The type is float.

ref_antenna_name

Required if cal_style = “sky”. Phase reference antenna.

sky_catalog

Required if cal_style = “sky”. Name of calibration catalog.

sky_field

Deprecated, will be removed in version 2.5. Only used if cal_style is ‘sky’. Short string describing field center or dominant source.

time_range

Time range (in JD) that cal solutions are valid for.list: [start_time, end_time] in JD. Should only be set if Ntimes is 1.

total_quality_array

Array of qualities of the calibration for entire arrays. The shape depends on cal_type, if the cal_type is ‘gain’, the shape is: (1, Nfreqs, Ntimes, Njones) or (Nfreqs, Ntimes, Njones) if future_array_shapes=True and wide_band=False, or (Nspws, Ntimes, Njones) if wide_band=True. If the cal_type is ‘delay’, the shape is (1, 1, Ntimes, Njones) or (Nspws, Ntimes, Njones) if future_array_shapes=True, type = float.

Methods

class pyuvdata.UVCal[source]

A class defining calibration solutions for interferometric data.

UVParameter objects

For full list see the documentation on ReadTheDocs: http://pyuvdata.readthedocs.io/en/latest/. Some are always required, some are required for certain cal_types and cal_styles and others are always optional.

static new(**kwargs)[source]

Create a new UVCal object.

All parameters are passed through to the new_uvcal() function.

Returns:

UVCal – A new UVCal object.

Other Parameters:
  • time_array (ndarray of float) – Array of times in Julian Date.

  • antenna_positions (ndarray of float or dict of ndarray of float) – Array of antenna positions in ECEF coordinates in meters. If a dict, keys can either be antenna numbers or antenna names, and values are position arrays. Keys are interpreted as antenna numbers if they are integers, otherwise they are interpreted as antenna names if strings. You cannot provide a mix of different types of keys.

  • telescope_location (ndarray of float or str) – Telescope location as an astropy EarthLocation object or MoonLocation object.

  • telescope_name (str) – Telescope name.

  • cal_style (str) – Calibration style. Options are ‘sky’ or ‘redundant’.

  • freq_array (ndarray of float, optional) – Array of frequencies in Hz. If given, the freq_range parameter is ignored, and the cal type is assumed to be “gain”.

  • freq_range (ndarray of float, optional) – Frequency ranges in Hz. Must be given if freq_array is not given. If given, the calibration is assumed to be wide band. The array shape should be (Nspws, 2)

  • gain_convention (str) – Gain convention. Options are ‘divide’ or ‘multiply’.

  • x_orientation (str) – Orientation of the x-axis. Options are ‘east’, ‘north’, ‘e’, ‘n’, ‘ew’, ‘ns’.

  • jones_array (ndarray of int or str) – Array of Jones polarization integers. If a string, options are ‘linear’ or ‘circular’ (which will be converted to the appropriate Jones integers as length-4 arrays).

  • cal_type (str, optional) – Calibration type. Options are ‘delay’, ‘gain’. Forced to be ‘gain’ if freq_array is given, and by default set to ‘delay’ if not.

  • integration_time (float or ndarray of float, optional) – Integration time in seconds. If not provided, it will be derived from the time_array, as the difference between successive times (with the last time-diff appended). If not provided and the number of unique times is one, then a warning will be raised and the integration time set to 1 second. If a float is provided, it will be used for all integrations. If an ndarray is provided, it must have the same shape as time_array.

  • channel_width (float or ndarray of float, optional) – Channel width in Hz. If not provided, it will be derived from the freq_array, as the difference between successive frequencies (with the last frequency-diff appended). If a float is provided, it will be used for all channels. If not provided and freq_array is length-one, the channel_width will be set to 1 Hz (and a warning issued). If an ndarray is provided, it must have the same shape as freq_array.

  • antenna_names (list of str, optional) – List of antenna names. If not provided, antenna numbers will be used to form the antenna_names, according to the antname_format. antenna_names need not be provided if antenna_positions is a dict with string keys.

  • antenna_numbers (list of int, optional) – List of antenna numbers. If not provided, antenna names will be used to form the antenna_numbers, but in this case the antenna_names must be strings that can be converted to integers. antenna_numbers need not be provided if antenna_positions is a dict with integer keys.

  • antname_format (str, optional) – Format string for antenna names. Default is ‘{0:03d}’.

  • ant_array (ndarray of int, optional) – Array of antenna numbers actually found in data (in the order of the data in gain_array etc.)

  • flex_spw_id_array (ndarray of int, optional) – Array of spectral window IDs. If not provided, it will be set to an array of zeros and only one spw will be used.

  • ref_antenna_name (str, optional) – Name of reference antenna. Only required for sky calibrations.

  • sky_catalog (str, optional) – Name of sky catalog. Only required for sky calibrations.

  • empty (bool, optional) – If True, create an empty UVCal object, i.e. add initialized data arrays (eg. gain_array). By default, this function creates a metadata-only object. You can pass in data arrays directly to the constructor if you want to create a fully-populated object.

  • data (dict of array_like, optional) – Dictionary containing optional data arrays. Possible keys are: ‘gain_array’, ‘delay_array’, ‘quality_array’, ‘flag_array’, and ‘total_quality_array’. If any entry is provided, the output will contain all necessary data-like arrays. Any key not provided in this case will be set to default “empty” values (e.g. all ones for gains, all False for flags) if they are required or None if they are not required (i.e. input_flag_array, quality_array, total_quality_array).

  • history (str, optional) – History string to be added to the object. Default is a simple string containing the date and time and pyuvdata version.

  • astrometry_library (str) – Library used for calculating LSTs. Allowed options are ‘erfa’ (which uses the pyERFA), ‘novas’ (which uses the python-novas library), and ‘astropy’ (which uses the astropy utilities). Default is erfa unless the telescope_location frame is MCMF (on the moon), in which case the default is astropy.

  • **kwargs – All other keyword arguments are added to the object as attributes.

property data_like_parameters

Iterate defined parameters which are data-like (not metadata-like).

property metadata_only

Property that determines whether this is a metadata only object.

An object is metadata only if data_array, nsample_array and flag_array are all None.

use_future_array_shapes()[source]

Change the array shapes of this object to match the planned future shapes.

This method sets allows users to convert to the planned array shapes changes before the changes go into effect. This method sets the future_array_shapes parameter on this object to True.

use_current_array_shapes()[source]

Change the array shapes of this object to match the current shapes.

This method sets allows users to convert back to the current array shapes. This method sets the future_array_shapes parameter on this object to False.

set_telescope_params(overwrite=False)[source]

Set telescope related parameters.

If the telescope_name is in the known_telescopes, set the telescope location to the value for the known telescope. Also set the antenna positions if they are not set on the object and are available for the telescope.

Parameters:

overwrite (bool) – Option to overwrite existing telescope-associated parameters with the values from the known telescope.

Raises:

ValueError – if the telescope_name is not in known telescopes

set_lsts_from_time_array(background=False, astrometry_library=None)[source]

Set the lst_array based from the time_array.

Parameters:

background (bool, False) – When set to True, start the calculation on a threading.Thread in the background and return the thread to the user.

Returns:

proc (None or threading.Thread instance) – When background is set to True, a thread is returned which must be joined before the lst_array exists on the UVCal object.

check(check_extra=True, run_check_acceptability=True, check_freq_spacing=False, lst_tol=3.63610260832152e-07)[source]

Add some extra checks on top of checks on UVBase class.

Check that required parameters exist. Check that parameters have appropriate shapes and optionally that the values are acceptable.

Parameters:
  • check_extra (bool) – If true, check all parameters, otherwise only check required parameters.

  • run_check_acceptability (bool) – Option to check if values in parameters are acceptable.

  • check_freq_spacing (bool) – Option to check if frequencies are evenly spaced and the spacing is equal to their channel_width. This is not required for UVCal objects in general but is required to write to calfits files.

  • lst_tol (float or None) – Tolerance level at which to test LSTs against their expected values. If provided as a float, must be in units of radians. If set to None, the default precision tolerance from the lst_array parameter is used (1 mas). Default value is 75 mas, which is set by the predictive uncertainty in IERS calculations of DUT1 (RMS is of order 1 ms, with with a 5-sigma threshold for detection is used to prevent false issues from being reported), which for some observatories sets the precision with which these values are written.

Returns:

bool – True if check passes

Raises:

ValueError – if parameter shapes or types are wrong or do not have acceptable values (if run_check_acceptability is True)

copy(metadata_only=False)[source]

Make and return a copy of the UVCal object.

Parameters:

metadata_only (bool) – If True, only copy the metadata of the object.

Returns:

UVCal – Copy of self.

ant2ind(antnum)[source]

Get the index in data arrays for an antenna number.

Parameters:

antnum (int) – Antenna number to get index for.

Returns:

int – Antenna index in data arrays.

jpol2ind(jpol)[source]

Get the index in data arrays for an antenna polarization.

Parameters:

jpol (int or str) – Antenna polarization to get index for.

Returns:

int – Antenna polarization index in data arrays

get_gains(ant, jpol=None, squeeze_pol=True)[source]

Get the gain associated with an antenna and/or polarization.

Parameters:
  • ant (int or length 2 tuple of ints or int and str) – Antenna or antenna and polarization to get gains for. If it’s a length 2 tuple, the second value must be an antenna polarization int or string parsable by jpol2ind.

  • jpol (int or str, optional) – Instrumental polarization to request. Ex. ‘Jxx’

  • squeeze_pol (bool) – Option to squeeze pol dimension if possible.

Returns:

complex ndarray – Gain solution of shape (Nfreqs, Ntimes, Njones) or (Nfreqs, Ntimes) if jpol is set or if squeeze_pol is True and Njones = 1.

get_flags(ant, jpol=None, squeeze_pol=True)[source]

Get the flags associated with an antenna and/or polarization.

Parameters:
  • ant (int or length 2 tuple of ints or int and str) – Antenna or antenna and polarization to get gains for. If it’s a length 2 tuple, the second value must be an antenna polarization int or string parsable by jpol2ind.

  • jpol (int or str, optional) – Instrumental polarization to request. Ex. ‘Jxx’

  • squeeze_pol (bool) – Option to squeeze pol dimension if possible.

Returns:

boolean ndarray – Flags of shape (Nfreqs, Ntimes, Njones) or (Nfreqs, Ntimes) if jpol is set or if squeeze_pol is True and Njones = 1.

get_quality(ant, jpol=None, squeeze_pol=True)[source]

Get the qualities associated with an antenna and/or polarization.

Parameters:
  • ant (int or length 2 tuple of ints or int and str) – Antenna or antenna and polarization to get gains for. If it’s a length 2 tuple, the second value must be an antenna polarization int or string parsable by jpol2ind.

  • jpol (int or str, optional) – Instrumental polarization to request. Ex. ‘Jxx’

  • squeeze_pol (bool) – Option to squeeze pol dimension if possible.

Returns:

float ndarray – Qualities of shape (Nfreqs, Ntimes, Njones) or (Nfreqs, Ntimes) if jpol is not None or if squeeze_pol is True and Njones = 1.

reorder_antennas(order='number', run_check=True, check_extra=True, run_check_acceptability=True)[source]

Arrange the antenna axis according to desired order.

Parameters:

order (str or array like of int) – If a string, allowed values are “name” and “number” to sort on the antenna name or number respectively. A ‘-’ can be prepended to signify descending order instead of the default ascending order (e.g. “-number”). An array of integers of length Nants representing indexes along the existing ant_array can also be supplied to sort in any desired order (note these are indices into the ant_array not antenna numbers).

Returns:

None

Raises:

ValueError – Raised if order is not an allowed string or is an array that does not contain all the required numbers.

reorder_freqs(spw_order=None, channel_order=None, select_spw=None, run_check=True, check_extra=True, run_check_acceptability=True)[source]

Arrange the frequency axis according to desired order.

Parameters:
  • spw_order (str or array_like of int) – A string describing the desired order of spectral windows along the frequecy axis. Allowed strings include number (sort on spectral window number) and freq (sort on median frequency). A ‘-’ can be prepended to signify descending order instead of the default ascending order, e.g., if you have SPW #1 and 2, and wanted them ordered as [2, 1], you would specify -number. Alternatively, one can supply an index array of length Nspws that specifies how to shuffle the spws (this is not the desired final spw order). Default is to apply no sorting of spectral windows.

  • channel_order (str or array_like of int) – A string describing the desired order of frequency channels within a spectral window. Allowed strings are “freq” and “-freq”, which will sort channels within a spectral window by ascending or descending frequency respectively. Alternatively, one can supply an index array of length Nfreqs that specifies the new order. Default is to apply no sorting of channels within a single spectral window. Note that proving an array_like of ints will cause the values given to spw_order and select_spw to be ignored.

  • select_spw (int or array_like of int) – An int or array_like of ints which specifies which spectral windows to apply sorting. Note that setting this argument will cause the value given to spw_order to be ignored.

Returns:

None

Raises:

ValueError – Raised if select_spw contains values not in spw_array, or if channel_order is not the same length as freq_array.

reorder_times(order='time', run_check=True, check_extra=True, run_check_acceptability=True)[source]

Arrange the time axis according to desired order.

Parameters:

order (str or array like of int) – If a string, allowed value is “time” or “-time” to sort on the time in ascending or descending order respectively. An array of integers of length Ntimes representing indexes along the existing time_array can also be supplied to sort in any desired order.

Returns:

None

Raises:

ValueError – Raised if order is not an allowed string or is an array that does not contain all the required indices.

reorder_jones(order='name', run_check=True, check_extra=True, run_check_acceptability=True)[source]

Arrange the jones element axis according to desired order.

Parameters:

order (str or array like of int) – If a string, allowed values are “name” and “number” to sort on the jones element name or number respectively. A ‘-’ can be prepended to signify descending order instead of the default ascending order (e.g. “-number”). An array of integers of length Njones representing indexes along the existing jones_array can also be supplied to sort in any desired order.

Returns:

None

Raises:

ValueError – Raised if order is not an allowed string or is an array that does not contain all the required indices.

convert_to_gain(freq_array=None, channel_width=None, delay_convention='minus', run_check=True, check_extra=True, run_check_acceptability=True)[source]

Convert non-gain cal_types to gains.

For the delay cal_type the gain is calculated as:

gain = 1 * exp((+/-) * 2 * pi * j * delay * frequency) where the (+/-) is dictated by the delay_convention

Parameters:
  • delay_convention (str) – Exponent sign to use in the conversion, can be “plus” or “minus”.

  • freq_array (array of float) – Frequencies to convert to gain at, units Hz. Not providing a freq_array is deprecated, but until version 3.0, if it is not provided and freq_array exists on the object, freq_array will be used.

  • run_check (bool) – Option to check for the existence and proper shapes of parameters after converting.

  • check_extra (bool) – Option to check optional parameters as well as required ones.

  • run_check_acceptability (bool) – Option to check acceptable range of the values of parameters after converting.

fast_concat(other, axis, inplace=False, verbose_history=False, run_check=True, check_extra=True, run_check_acceptability=True)[source]

Concatenate two UVCal objects along specified axis with almost no checking.

Warning! This method assumes all the metadata along other axes is sorted the same way. The __add__ method is much safer, it checks all the metadata, but it is slower. Some quick checks are run, but this method doesn’t make any guarantees that the resulting object is correct.

Parameters:
  • other (UVCal object or list of UVCal objects) – UVCal object or list of UVCal objects which will be added to self.

  • axis (str) – Axis to concatenate files along. This enables fast concatenation along the specified axis without the normal checking that all other metadata agrees. Allowed values are: ‘antenna’, ‘time’, ‘freq’, ‘spw’, ‘jones’ (‘freq’ is not allowed for delay or wideband objects and ‘spw’ is only allowed for wideband objects).

  • inplace (bool) – If True, overwrite self as we go, otherwise create a third object as the sum of the two.

  • verbose_history (bool) – Option to allow more verbose history. If True and if the histories for the two objects are different, the combined object will keep all the history of both input objects (if many objects are combined in succession this can lead to very long histories). If False and if the histories for the two objects are different, the combined object will have the history of the first object and only the parts of the second object history that are unique (this is done word by word and can result in hard to interpret histories).

  • run_check (bool) – Option to check for the existence and proper shapes of parameters after combining objects.

  • check_extra (bool) – Option to check optional parameters as well as required ones.

  • run_check_acceptability (bool) – Option to check acceptable range of the values of parameters after combining objects.

select(antenna_nums=None, antenna_names=None, frequencies=None, freq_chans=None, spws=None, times=None, jones=None, run_check=True, check_extra=True, run_check_acceptability=True, inplace=True)[source]

Downselect data to keep on the object along various axes.

Axes that can be selected along include antennas, frequencies, times and antenna polarization (jones).

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

Parameters:
  • antenna_nums (array_like of int, optional) – The antennas numbers to keep in the object (antenna positions and names for the removed antennas will be retained). This cannot be provided if antenna_names is also provided.

  • antenna_names (array_like of str, optional) – The antennas names to keep in the object (antenna positions and names for the removed antennas will be retained). This cannot be provided if antenna_nums is also provided.

  • frequencies (array_like of float, optional) – The frequencies to keep in the object, each value passed here should exist in the freq_array.

  • freq_chans (array_like of int, optional) – The frequency channel numbers to keep in the object.

  • spws (array_like of in, optional) – The spectral window numbers to keep in the object. If this is not a wide-band object and frequencies or freq_chans is not None, frequencies that match any of the specifications will be kept (i.e. the selections will be OR’ed together).

  • times (array_like of float, optional) – The times to keep in the object, each value passed here should exist in the time_array.

  • jones (array_like of int or str, optional) – The antenna polarizations numbers to keep in the object, each value passed here should exist in the jones_array. If passing strings, the canonical polarization strings (e.g. “Jxx”, “Jrr”) are supported and if the x_orientation attribute is set, the physical dipole strings (e.g. “Jnn”, “Jee”) are also supported.

  • run_check (bool) – Option to check for the existence and proper shapes of parameters after downselecting data on this object (the default is True, meaning the check will be run).

  • check_extra (bool) – Option to check optional parameters as well as required ones (the default is True, meaning the optional parameters will be checked).

  • run_check_acceptability (bool) – Option to check acceptable range of the values of parameters after downselecting data on this object (the default is True, meaning the acceptable range check will be done).

  • inplace (bool) – Option to perform the select directly on self or return a new UVCal object with just the selected data (the default is True, meaning the select will be done on self).

classmethod initialize_from_uvdata(uvdata, gain_convention, cal_style, future_array_shapes=True, metadata_only=True, times=None, frequencies=None, jones=None, **kwargs)[source]

Initialize this object based on a UVData object.

Internally, this function takes whatever default parameters it can from the UVData object, and then overwrites them with any parameters passed in as kwargs. It then uses new_uvcal() to construct the UVCal object.

Parameters:
  • uvdata (UVData object) – The UVData object to initialize from.

  • gain_convention (str) – Gain convention. Options are ‘divide’ or ‘multiply’.

  • cal_style (str) – Calibration style. Options are ‘sky’ or ‘redundant’.

  • future_array_shapes (bool) – Option to use the future array shapes (see use_future_array_shapes for details). Note that this option is deprecated and will be removed in v3.

  • metadata_only (bool) – Option to only initialize the metadata. If False, this method also initializes the data-like arrays to zeros/ones as appropriate (or False for the flag_array) with the appropriate sizes.

  • times (array_like of float, optional) – Deprecated alias for time_array. Will be removed in v2.5.

  • frequencies (array_like of float, optional) – Deprecated alias for freq_array. Will be removed in v2.5.

  • jones (array_like of int, optional) – Deprecated alias for jones_array. Will be removed in v2.5.

Returns:

UVCal – A new UVCal object with default parameters.

Other Parameters:
  • cal_type (str) – Calibration type. Must be ‘gain’ or ‘delay’. Default ‘gain’. Note that in contrast to new_uvcal(), this function sets the default to ‘gain’, instead of trying to set it based on the input freq_range.

  • jones_array (array_like of int or str, optional) – Array of Jones polarizations. Taken from UVData object is possible (i.e. it is single-polarization data), otherwise must be specified. May be specified as ‘circular’ or ‘linear’ to indicate all circular or all linear polarizations.

  • wide_band (bool) – Whether the calibration is wide-band (i.e. one calibration parameter per spectral window instead of per channel). Automatically set to True if cal_type is ‘delay’, otherwise the default is False. Note that if wide_band is True, the freq_range parameter is required.

  • include_uvdata_history (bool) – Whether to include the UVData object’s history in the UVCal object’s history. Default is True.

  • spw_array (array_like of int, optional) – Array of spectral window numbers. This will be taken from the UVData object if at all possible. The only time it is useful to pass this in explicitly is if wide_band=True and the spectral windows desired are not just (0, ..., Nspw-1).

  • time_array (ndarray of float) – Array of times in Julian Date.

  • antenna_positions (ndarray of float or dict of ndarray of float) – Array of antenna positions in ECEF coordinates in meters. If a dict, keys can either be antenna numbers or antenna names, and values are position arrays. Keys are interpreted as antenna numbers if they are integers, otherwise they are interpreted as antenna names if strings. You cannot provide a mix of different types of keys.

  • telescope_location (ndarray of float or str) – Telescope location as an astropy EarthLocation object or MoonLocation object.

  • telescope_name (str) – Telescope name.

  • freq_array (ndarray of float, optional) – Array of frequencies in Hz. If given, the freq_range parameter is ignored, and the cal type is assumed to be “gain”.

  • freq_range (ndarray of float, optional) – Frequency ranges in Hz. Must be given if freq_array is not given. If given, the calibration is assumed to be wide band. The array shape should be (Nspws, 2)

  • x_orientation (str) – Orientation of the x-axis. Options are ‘east’, ‘north’, ‘e’, ‘n’, ‘ew’, ‘ns’.

  • integration_time (float or ndarray of float, optional) – Integration time in seconds. If not provided, it will be derived from the time_array, as the difference between successive times (with the last time-diff appended). If not provided and the number of unique times is one, then a warning will be raised and the integration time set to 1 second. If a float is provided, it will be used for all integrations. If an ndarray is provided, it must have the same shape as time_array.

  • channel_width (float or ndarray of float, optional) – Channel width in Hz. If not provided, it will be derived from the freq_array, as the difference between successive frequencies (with the last frequency-diff appended). If a float is provided, it will be used for all channels. If not provided and freq_array is length-one, the channel_width will be set to 1 Hz (and a warning issued). If an ndarray is provided, it must have the same shape as freq_array.

  • antenna_names (list of str, optional) – List of antenna names. If not provided, antenna numbers will be used to form the antenna_names, according to the antname_format. antenna_names need not be provided if antenna_positions is a dict with string keys.

  • antenna_numbers (list of int, optional) – List of antenna numbers. If not provided, antenna names will be used to form the antenna_numbers, but in this case the antenna_names must be strings that can be converted to integers. antenna_numbers need not be provided if antenna_positions is a dict with integer keys.

  • antname_format (str, optional) – Format string for antenna names. Default is ‘{0:03d}’.

  • ant_array (ndarray of int, optional) – Array of antenna numbers actually found in data (in the order of the data in gain_array etc.)

  • flex_spw_id_array (ndarray of int, optional) – Array of spectral window IDs. If not provided, it will be set to an array of zeros and only one spw will be used.

  • ref_antenna_name (str, optional) – Name of reference antenna. Only required for sky calibrations.

  • sky_catalog (str, optional) – Name of sky catalog. Only required for sky calibrations.

  • empty (bool, optional) – If True, create an empty UVCal object, i.e. add initialized data arrays (eg. gain_array). By default, this function creates a metadata-only object. You can pass in data arrays directly to the constructor if you want to create a fully-populated object.

  • data (dict of array_like, optional) – Dictionary containing optional data arrays. Possible keys are: ‘gain_array’, ‘delay_array’, ‘quality_array’, ‘flag_array’, and ‘total_quality_array’. If any entry is provided, the output will contain all necessary data-like arrays. Any key not provided in this case will be set to default “empty” values (e.g. all ones for gains, all False for flags) if they are required or None if they are not required (i.e. input_flag_array, quality_array, total_quality_array).

  • history (str, optional) – History string to be added to the object. Default is a simple string containing the date and time and pyuvdata version.

  • astrometry_library (str) – Library used for calculating LSTs. Allowed options are ‘erfa’ (which uses the pyERFA), ‘novas’ (which uses the python-novas library), and ‘astropy’ (which uses the astropy utilities). Default is erfa unless the telescope_location frame is MCMF (on the moon), in which case the default is astropy.

  • **kwargs – All other keyword arguments are added to the object as attributes.

read_calfits(filename, **kwargs)[source]

Read in data from calfits file(s).

Parameters:
  • filename (str) – The calfits file to read from.

  • read_data (bool) – Read in the gains or delays, quality arrays and flag arrays. If set to False, only the metadata will be read in. Setting read_data to False results in a metadata only object.

  • background_lsts (bool) – When set to True, the lst_array is calculated in a background thread.

  • run_check (bool) – Option to check for the existence and proper shapes of parameters after reading in the file.

  • check_extra (bool) – Option to check optional parameters as well as required ones.

  • run_check_acceptability (bool) – Option to check acceptable range of the values of parameters after reading in the file.

  • use_future_array_shapes (bool) – Option to convert to the future planned array shapes before the changes go into effect by removing the spectral window axis.

  • astrometry_library (str) – Library used for calculating LSTs. Allowed options are ‘erfa’ (which uses the pyERFA), ‘novas’ (which uses the python-novas library), and ‘astropy’ (which uses the astropy utilities). Default is erfa unless the telescope_location frame is MCMF (on the moon), in which case the default is astropy.

read_fhd_cal(cal_file, obs_file, layout_file=None, settings_file=None, **kwargs)[source]

Read data from an FHD cal.sav file.

Parameters:
  • cal_file (str or list of str) – The cal.sav file or list of files to read from.

  • obs_file (str or list of str) – The obs.sav file or list of files to read from.

  • layout_file (str) – The FHD layout file. Required for antenna_positions to be set.

  • settings_file (str or list of str, optional) – The settings_file or list of files to read from. Optional, but very useful for provenance.

  • raw (bool) – Option to use the raw (per antenna, per frequency) solution or to use the fitted (polynomial over phase/amplitude) solution. Default is True (meaning use the raw solutions).

  • read_data (bool) – Read in the gains, quality array and flag data. If set to False, only the metadata will be read in. Setting read_data to False results in a metadata only object. If read_data is False, a settings file must be provided.

  • background_lsts (bool) – When set to True, the lst_array is calculated in a background thread.

  • extra_history (str or list of str, optional) – String(s) to add to the object’s history parameter.

  • run_check (bool) – Option to check for the existence and proper shapes of parameters after reading in the file.

  • check_extra (bool) – Option to check optional parameters as well as required ones.

  • run_check_acceptability (bool) – Option to check acceptable range of the values of parameters after reading in the file.

  • use_future_array_shapes (bool) – Option to convert to the future planned array shapes before the changes go into effect by removing the spectral window axis.

  • astrometry_library (str) – Library used for calculating LSTs. Allowed options are ‘erfa’ (which uses the pyERFA), ‘novas’ (which uses the python-novas library), and ‘astropy’ (which uses the astropy utilities). Default is erfa unless the telescope_location frame is MCMF (on the moon), in which case the default is astropy.

read(filename, *, axis=None, file_type=None, read_data=True, use_future_array_shapes=False, run_check=True, check_extra=True, run_check_acceptability=True, obs_file=None, layout_file=None, settings_file=None, raw=True, extra_history=None)[source]

Read a generic file into a UVCal object.

This method supports a number of different types of files. Universal parameters (required and optional) are listed directly below, followed by parameters used by all file types related to checking. Each file type also has its own set of optional parameters that are listed at the end of this docstring.

Parameters:
  • filename (str or array_like of str) – The file(s) or list(s) (or array(s)) of files to read from.

  • file_type (str) – One of [‘calfits’, ‘fhd’] or None. If None, the code attempts to guess what the file type is based on file extensions (FHD: .sav, .txt; uvfits: .calfits). Note that if a list of datasets is passed, the file type is determined from the first dataset.

  • axis (str) – Axis to concatenate files along. This enables fast concatenation along the specified axis without the normal checking that all other metadata agrees. This method does not guarantee correct resulting objects. Please see the docstring for fast_concat for details. Allowed values are: ‘antenna’, ‘time’, ‘freq’, ‘spw’, ‘jones’ (‘freq’ is not allowed for delay or wideband objects and ‘spw’ is only allowed for wideband objects). Only used if multiple files are passed.

  • read_data (bool) – Read in the gains or delays, quality arrays and flag arrays. If set to False, only the metadata will be read in. Setting read_data to False results in a metadata only object.

  • use_future_array_shapes (bool) – Option to convert to the future planned array shapes before the changes go into effect by removing the spectral window axis.

Checking:
  • run_check (bool) – Option to check for the existence and proper shapes of parameters after reading in the file.

  • check_extra (bool) – Option to check optional parameters as well as required ones.

  • run_check_acceptability (bool) – Option to check acceptable range of the values of parameters after reading in the file.

FHD:
  • obs_file (str or list of str) – The obs.sav file or list of files to read from. This is required for FHD files.

  • layout_file (str) – The FHD layout file. Required for antenna_positions to be set.

  • settings_file (str or list of str, optional) – The settings_file or list of files to read from. Optional, but very useful for provenance.

  • raw (bool) – Option to use the raw (per antenna, per frequency) solution or to use the fitted (polynomial over phase/amplitude) solution. Default is True (meaning use the raw solutions).

classmethod from_file(filename, *, axis=None, file_type=None, read_data=True, use_future_array_shapes=False, run_check=True, check_extra=True, run_check_acceptability=True, obs_file=None, layout_file=None, settings_file=None, raw=True, extra_history=None)[source]

Initialize a new UVCal object by reading the input file.

This method supports a number of different types of files. Universal parameters (required and optional) are listed directly below, followed by parameters used by all file types related to checking. Each file type also has its own set of optional parameters that are listed at the end of this docstring.

Parameters:
  • filename (str or array_like of str) – The file(s) or list(s) (or array(s)) of files to read from.

  • file_type (str) – One of [‘calfits’, ‘fhd’] or None. If None, the code attempts to guess what the file type is based on file extensions (FHD: .sav, .txt; uvfits: .calfits). Note that if a list of datasets is passed, the file type is determined from the first dataset.

  • axis (str) – Axis to concatenate files along. This enables fast concatenation along the specified axis without the normal checking that all other metadata agrees. This method does not guarantee correct resulting objects. Please see the docstring for fast_concat for details. Allowed values are: ‘antenna’, ‘time’, ‘freq’, ‘spw’, ‘jones’ (‘freq’ is not allowed for delay or wideband objects and ‘spw’ is only allowed for wideband objects). Only used if multiple files are passed.

  • read_data (bool) – Read in the gains or delays, quality arrays and flag arrays. If set to False, only the metadata will be read in. Setting read_data to False results in a metadata only object.

  • use_future_array_shapes (bool) – Option to convert to the future planned array shapes before the changes go into effect by removing the spectral window axis.

Checking:
  • run_check (bool) – Option to check for the existence and proper shapes of parameters after reading in the file.

  • check_extra (bool) – Option to check optional parameters as well as required ones.

  • run_check_acceptability (bool) – Option to check acceptable range of the values of parameters after reading in the file.

FHD:
  • obs_file (str or list of str) – The obs.sav file or list of files to read from. This is required for FHD files.

  • layout_file (str) – The FHD layout file. Required for antenna_positions to be set.

  • settings_file (str or list of str, optional) – The settings_file or list of files to read from. Optional, but very useful for provenance.

  • raw (bool) – Option to use the raw (per antenna, per frequency) solution or to use the fitted (polynomial over phase/amplitude) solution. Default is True (meaning use the raw solutions).

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

Write the data to a calfits file.

Parameters:
  • filename (str) – The calfits file to write to.

  • run_check (bool) – Option to check for the existence and proper shapes of parameters before writing the file.

  • check_extra (bool) – Option to check optional parameters as well as required ones.

  • run_check_acceptability (bool) – Option to check acceptable range of the values of parameters before writing the file.

  • clobber (bool) – Option to overwrite the filename if the file already exists.

Raises:

ValueError – If the UVCal object is a metadata only object.

last updated: 2024-03-25