Full class list and documentation

Primary Beams

Simulations may be run using pyuvdata UVBeam objects or AnalyticBeam objects.

class pyuvsim.AnalyticBeam(type, sigma=None, diameter=None)[source]

Defines an object with similar functionality to pyuvdata.UVBeam

Directly calculates jones matrices at given azimuths and zenith angles from analytic functions.

Supports uniform (unit response in all directions), gaussian, and Airy function beam types.

interp(az_array, za_array, freq_array)[source]

Evaluate the primary beam at given az, za locations (in radians).

(similar to UVBeam.interp)

Parameters:
  • az_array – az values to evaluate at in radians (same length as za_array) The azimuth here has the UVBeam convention: North of East(East=0, North=pi/2)
  • za_array – za values to evaluate at in radians (same length as az_array)
  • freq_array – frequency values to evaluate at
Returns:

an array of beam 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)

Antenna objects

class pyuvsim.Antenna(name, number, enu_position, beam_id)[source]

Defines an object that can return a Jones matrix and has a specified location, name, and number.

One of these defined for each antenna in the array.

get_beam_jones(array, source_alt_az, frequency)[source]

2x2 array of Efield vectors in Az/Alt

Parameters:
  • array – az values to evaluate (same length as za_array)
  • source_az_za – Tuple or list (azimuth, zenith angle) in radians
  • frequency – (float) frequency in Hz
Returns:

jones_matrix, A 2x2 float array. The first axis is feed, the

second axis is vector component on the sky in az/za.

Profiling

Developers will want to use line profiling tools to track down bottlenecks and find room for performance improvements. These tools better allow profilers to operate within the MPI environment.

Use the line profiler when requested.

If the set_profiler method is called, the profile decorators throughout the code will do line profiling on all instance methods.

pyuvsim.profiling.set_profiler()[source]

If profiling is requested, then assign it to the builtins

Simulation setup

Simulations are run from .yaml files with specified keywords. See usage.

pyuvsim.simsetup.create_mock_catalog(time, arrangement='zenith', array_location=None, Nsrcs=None, alt=None, save=False, min_alt=None, rseed=None)[source]

Create a mock catalog.

Sources are defined in an AltAz frame at the given time, then returned in ICRS ra/dec coordinates.

Parameters:
  • time (float or astropy Time object) – Julian date
  • arrangement (str) –

    Point source pattern (default = 1 source at zenith). Accepted arrangements:

    triangle: Three point sources forming a triangle around the zenith
    cross: An asymmetric cross
    zenith: Some number of sources placed at the zenith.
    off-zenith: A single source off zenith
    long-line: Horizon to horizon line of point sources
    hera_text: Spell out HERA around the zenith
    random: Randomly distributed point sources near zenith
  • Nsrcs (int) – Number of sources to put at zenith
  • array_location (EarthLocation object) – [Default = HERA site]
  • alt (float) – For off-zenith and triangle arrangements, altitude to place sources. (deg)
  • save (bool) – Save mock catalog as npz file.
  • rseed (int) – If using the random configuration, pass in a RandomState seed.
Returns:

List of pyuvsim.Source objects mock_kwds: (dictionary) The keywords defining this source catalog

Return type:

catalog

pyuvsim.simsetup.initialize_uvdata_from_params(obs_params)[source]

Construct a uvdata object from parameters in a valid yaml file.

Sufficient information must be provided by the parameters to define time and frequency arrays and verify the channel widths and time steps. This will error if insufficient or incompatible parameters are defined.

The parameter dictionary may contain any valid UVData attributes as well.

Parameters:obs_params – Either an obs_param file name or a dictionary of parameters read in. Any uvdata parameters may be passed in through here.
Returns:uv_obj, beam_list, beam_dict, beam_ids
pyuvsim.simsetup.read_text_catalog(catalog_csv)[source]

Read in a text file of sources.

Args: catalog_csv: csv file with the following expected columns:

Source_ID: source id as a string of maximum 10 characters
ra_j2000: right ascension at J2000 epoch, in decimal degrees
dec_j2000: declination at J2000 epoch, in decimal degrees
flux_density_I: Stokes I flux density in Janskys
frequency: reference frequency (for future spectral indexing) [Hz] For now, all sources are flat spectrum.
Returns:List of pyuvsim.Source objects
pyuvsim.simsetup.read_votable_catalog(gleam_votable)[source]

Creates a list of pyuvsim source objects from a votable catalog.

Parameters:gleam_votable – Path to votable catalog file.
Returns:List of pyuvsim.Source objects

Tested on: GLEAM EGC catalog, version 2

pyuvsim.simsetup.uvdata_to_config_file(uvdata_in, param_filename=None, telescope_config_name='', layout_csv_name='', catalog='mock', path_out='.')[source]

Extract simulation configuration settings from uvfits.

Parameters:uvdata_in (UVData) – uvdata object.
Keywords:
param_filename (str, optional): output param file name, defaults to obsparam_#.yaml. telescope_config_name (str, optional): Name of yaml file file. Defaults to blank string. layout_csv_name (str, optional): Name of layout csv file. Defaults to blank string. catalog (str, optional): Path to catalog file, defaults to ‘mock’. path_out (str, optional): Where to put config files.
pyuvsim.simsetup.uvdata_to_telescope_config(uvdata_in, beam_filepath, layout_csv_name=None, telescope_config_name=None, return_names=False, path_out='.')[source]

From a uvfits file, generate telescope parameters files.

Output config files are written to the current directory, unless keep_path is set.

Parameters:
  • uvdata_in (UVData) – object to process
  • path_out (str) – Target directory for the config file.
  • beam_filepath (str) – Path to a beamfits file.
  • layout_csv_name (str, optional) – The name for the antenna positions csv file (Default <telescope_name>_layout.csv)
  • telescope_config_name (str, optional) – The name for the telescope config file (Default teleconfig_#number.yaml)
  • return_names (bool, optional) – Return the file names for loopback tests.
Returns:

if return_names, returns tuple (path, telescope_config_name, layout_csv_name)

pyuvsim.simsetup.write_catalog_to_file(filename, catalog)[source]

Writes out a catalog to a text file, readable with simsetup.read_catalog_text()

Parameters:
  • filename – Path to output file (string)
  • catalog – List of pyuvsim.Source objects
pyuvsim.simsetup.write_uvfits(uv_obj, param_dict, return_filename=False, dryrun=False)[source]

Parse output file information from parameters and write uvfits to file.

Parameters:
  • uv_obj – UVData object to write out.
  • param_dict – parameter dictionary defining output path, filename, and whether or not to clobber.
  • return_filename – (Default false) Return the file path
  • dryrun – (Default false) Don’t write to file.
Returns:

File path, if return_filename is True

Sources

class pyuvsim.Source(name, ra, dec, freq, stokes, pos_tol=2.220446049250313e-16)[source]

Defines a single point source at a given ICRS ra/dec coordinate, with a flux density defined by stokes parameters.

Used to calculate local coherency matrix for source brightness in AltAz frame at a specified time.

alt_az_calc(time, telescope_location)[source]

calculate the altitude & azimuth for this source at a time & location

Parameters:
  • time – astropy Time object
  • telescope_location – astropy EarthLocation object
Returns:

(altitude, azimuth) in radians

coherency_calc(time, telescope_location)[source]

Calculate the local coherency in alt/az basis for this source at a time & location.

The coherency is a 2x2 matrix giving electric field correlation in Jy. It’s specified on the object as a coherency in the ra/dec basis, but must be rotated into local alt/az.

Parameters:
  • time – astropy Time object
  • telescope_location – astropy EarthLocation object
Returns:

local coherency in alt/az basis

pos_lmn(time, telescope_location)[source]

calculate the direction cosines of this source at a time & location

Parameters:
  • time – astropy Time object
  • telescope_location – astropy EarthLocation object
Returns:

(l, m, n) direction cosine values

Telescopes

Shared properties of all antennas, such as array center location and primary beam objects.

class pyuvsim.Telescope(telescope_name, telescope_location, beam_list)[source]

Defines the location and name of the observing site, and holds the list of beam objects used by the array.

UV Simulation functions

Methods for setting up the simulation, handling concurrent MPI processes, and saving computed visibilities.

pyuvsim.uvsim.uvdata_to_task_list(input_uv, sources, beam_list, beam_dict=None)[source]

Create task list from pyuvdata compatible input file.

Parameters:
  • input_uv (UVData) – UVData object to use
  • sources – array of Source objects
  • beam_list – (list of UVBeam or AnalyticBeam objects
  • beam_dict (dict, optional) – dict mapping antenna number to beam index in beam_list
Returns:

List of task parameters to be send to UVEngines with the task parameters defined in UVTask objectself. This function extracts time, freq, Antenna1, Antenna2.

pyuvsim.uvsim.run_uvsim(input_uv, beam_list, beam_dict=None, catalog_file=None, mock_keywords=None, uvdata_file=None, obs_param_file=None, telescope_config_file=None, antenna_location_file=None)[source]

Run uvsim

Parameters:
  • input_uv – An input UVData object, containing baseline/time/frequency information.
  • beam_list – A list of UVBeam and/or AnalyticBeam objects
  • beam_dict – Dictionary of {antenna_name : beam_ID}, where beam_id is an index in the beam_list. This assigns beams to antennas. Default: All antennas get the 0th beam in the beam_list.
  • catalog_file – Catalog file name. Default: Create a mock catalog
  • mock_keywords – Settings for a mock catalog (see keywords of create_mock_catalog)
  • uvdata_file – Name of input UVData file if running from a file.
  • obs_param_file – Parameter filename if running from config files.
  • telescope_config_file – Telescope configuration file if running from config files.
  • antenna_location_file – antenna_location file if running from config files.
pyuvsim.uvsim.initialize_uvdata(uvtask_list, source_list_name, uvdata_file=None, obs_param_file=None, telescope_config_file=None, antenna_location_file=None)[source]

Initialize an empty uvdata object to fill with simulation.

Parameters:
  • uvtask_list – List of uvtasks to simulate.
  • source_list_name – Name of source list file or mock catalog.
  • uvdata_file – Name of input UVData file or None if initializing from config files.
  • obs_param_file – Name of observation parameter config file or None if initializing from a UVData file.
  • telescope_config_file – Name of telescope config file or None if initializing from a UVData file.
  • antenna_location_file – Name of antenna location file or None if initializing from a UVData file.
pyuvsim.uvsim.serial_gather(uvtask_list, uv_out)[source]

Loop over uvtask list, acquire visibilities and add to uvdata object.