Radio Interferometer Measurement Equation¶
Functions used to compute the terms of the Radio Interferometer Measurement Equation (RIME). It describes the response of an interferometer to a sky model.
where for antenna \(p\) and \(q\), and source \(s\):
- \(G_{p}\) represents direction-independent effects.
- \(E_{ps}\) represents direction-dependent effects.
- \(L_{p}\) represents the feed rotation.
- \(K_{ps}\) represents the phase delay term.
- \(B_{s}\) represents the brightness matrix.
The RIME is more formally described in the following four papers:
- I. A full-sky Jones formalism
- II. Calibration and direction-dependent effects
- III. Addressing direction-dependent effects in 21cm WSRT observations of 3C147
- IV. A generalized tensor formalism
Numpy¶
predict_vis (time_index, antenna1, antenna2) |
Multiply Jones terms together to form model visibilities according to the following formula: |
phase_delay (lm, uvw, frequency[, convention]) |
Computes the phase delay (K) term: |
parallactic_angles (times, antenna_positions, …) |
Computes parallactic angles per timestep for the given reference antenna position and field centre. |
feed_rotation (parallactic_angles[, feed_type]) |
Computes the 2x2 feed rotation (L) matrix from the parallactic_angles . |
transform_sources (lm, parallactic_angles, …) |
Creates beam sampling coordinates suitable for use in beam_cube_dde() by: |
beam_cube_dde (beam, beam_lm_extents, …) |
Evaluates Direction Dependent Effects along a source’s path by interpolating the values of a complex beam cube at the source location. |
zernike_dde (coords, coeffs, noll_index) |
Computes Direction Dependent Effects by evaluating Zernicke Polynomials defined by coefficients coeffs and noll indexes noll_index at the specified coordinates coords . |
wsclean_predict (uvw, lm, source_type, flux, …) |
Predict visibilities from a WSClean sky model. |
-
africanus.rime.
predict_vis
(time_index, antenna1, antenna2, dde1_jones=None, source_coh=None, dde2_jones=None, die1_jones=None, base_vis=None, die2_jones=None)[source]¶ Multiply Jones terms together to form model visibilities according to the following formula:
\[V_{pq} = G_{p} \left( B_{pq} + \sum_{s} E_{ps} X_{pqs} E_{qs}^H \right) G_{q}^H\]where for antenna \(p\) and \(q\), and source \(s\):
- \(B_{{pq}}\) represent base coherencies.
- \(E_{{ps}}\) represents Direction-Dependent Jones terms.
- \(X_{{pqs}}\) represents a coherency matrix (per-source).
- \(G_{{p}}\) represents Direction-Independent Jones terms.
Generally, \(E_{ps}\), \(G_{p}\), \(X_{pqs}\) should be formed by using the RIME API functions and combining them together with
einsum()
.Please read the Notes
Parameters: time_index :
numpy.ndarray
Time index used to look up the antenna Jones index for a particular baseline with shape
(row,)
. Obtainable vianp.unique(time, return_inverse=True)[1]
.antenna1 :
numpy.ndarray
Antenna 1 index used to look up the antenna Jones for a particular baseline. with shape
(row,)
.antenna2 :
numpy.ndarray
Antenna 2 index used to look up the antenna Jones for a particular baseline. with shape
(row,)
.dde1_jones :
numpy.ndarray
, optional\(E_{ps}\) Direction-Dependent Jones terms for the first antenna. shape
(source,time,ant,chan,corr_1,corr_2)
source_coh :
numpy.ndarray
, optional\(X_{pqs}\) Direction-Dependent Coherency matrix for the baseline. with shape
(source,row,chan,corr_1,corr_2)
dde2_jones :
numpy.ndarray
, optional\(E_{qs}\) Direction-Dependent Jones terms for the second antenna. This is usually the same array as
dde1_jones
as this preserves the symmetry of the RIME.predict_vis
will perform the conjugate transpose internally. shape(source,time,ant,chan,corr_1,corr_2)
die1_jones :
numpy.ndarray
, optional\(G_{ps}\) Direction-Independent Jones terms for the first antenna of the baseline. with shape
(time,ant,chan,corr_1,corr_2)
base_vis :
numpy.ndarray
, optional\(B_{pq}\) base coherencies, added to source coherency summation before multiplication with die1_jones and die2_jones. shape
(row,chan,corr_1,corr_2)
.die2_jones :
numpy.ndarray
, optional\(G_{ps}\) Direction-Independent Jones terms for the second antenna of the baseline. This is usually the same array as
die1_jones
as this preserves the symmetry of the RIME.predict_vis
will perform the conjugate transpose internally. shape(time,ant,chan,corr_1,corr_2)
Returns: visibilities :
numpy.ndarray
Model visibilities of shape
(row,chan,corr_1,corr_2)
Notes
- Direction-Dependent terms (dde{1,2}_jones) and Independent (die{1,2}_jones) are optional, but if one is present, the other must be present.
- The inputs to this function involve
row
,time
andant
(antenna) dimensions. - Each
row
is associated with a pair of antenna Jones matrices at a particular timestep via thetime_index
,antenna1
andantenna2
inputs. - The
row
dimension must be an increasing partial order in time.
-
africanus.rime.
phase_delay
(lm, uvw, frequency, convention='fourier')[source]¶ Computes the phase delay (K) term:
\[ \begin{align}\begin{aligned}& {\Large e^{-2 \pi i (u l + v m + w (n - 1))} }\\& \textrm{where } n = \sqrt{1 - l^2 - m^2}\end{aligned}\end{align} \]Parameters: lm :
numpy.ndarray
LM coordinates of shape
(source, 2)
with L and M components in the last dimension.uvw :
numpy.ndarray
UVW coordinates of shape
(row, 3)
with U, V and W components in the last dimension.frequency :
numpy.ndarray
frequencies of shape
(chan,)
convention : {‘fourier’, ‘casa’}
Uses the \(e^{-2 \pi \mathit{i}}\) sign convention if
fourier
and \(e^{2 \pi \mathit{i}}\) ifcasa
.Returns: complex_phase :
numpy.ndarray
complex of shape
(source, row, chan)
Notes
Corresponds to the complex exponential of the Van Cittert-Zernike Theorem.
MeqTrees uses the CASA sign convention.
-
africanus.rime.
parallactic_angles
(times, antenna_positions, field_centre, backend='casa')[source]¶ Computes parallactic angles per timestep for the given reference antenna position and field centre.
Parameters: times :
numpy.ndarray
Array of Mean Julian Date times in seconds with shape
(time,)
,antenna_positions :
numpy.ndarray
Antenna positions of shape
(ant, 3)
in metres in the ITRF frame.field_centre :
numpy.ndarray
Field centre of shape
(2,)
in radiansbackend : {‘casa’, ‘test’}, optional
Backend to use for calculating the parallactic angles.
casa
defers to an implementation depending onpython-casacore
. This backend should be used by default.test
creates parallactic angles by multiplying thetimes
andantenna_position
arrays. It exist solely for testing.
Returns: parallactic_angles :
numpy.ndarray
Parallactic angles of shape
(time,ant)
-
africanus.rime.
feed_rotation
(parallactic_angles, feed_type='linear')[source]¶ Computes the 2x2 feed rotation (L) matrix from the
parallactic_angles
.\[\begin{split}\textrm{linear} \begin{bmatrix} cos(pa) & sin(pa) \\ -sin(pa) & cos(pa) \end{bmatrix} \qquad \textrm{circular} \begin{bmatrix} e^{-i pa} & 0 \\ 0 & e^{i pa} \end{bmatrix}\end{split}\]Parameters: parallactic_angles :
numpy.ndarray
floating point parallactic angles. Of shape
(pa0, pa1, ..., pan)
.feed_type : {‘linear’, ‘circular’}
The type of feed
Returns: feed_matrix :
numpy.ndarray
Feed rotation matrix of shape
(pa0, pa1,...,pan,2,2)
-
africanus.rime.
transform_sources
(lm, parallactic_angles, pointing_errors, antenna_scaling, frequency, dtype=None)[source]¶ Creates beam sampling coordinates suitable for use in
beam_cube_dde()
by:- Rotating
lm
coordinates by theparallactic_angles
- Adding
pointing_errors
- Scaling by
antenna_scaling
Parameters: lm :
numpy.ndarray
LM coordinates of shape
(src,2)
in radians offset from the phase centre.parallactic_angles :
numpy.ndarray
parallactic angles of shape
(time, antenna)
in radians.pointing_errors :
numpy.ndarray
LM pointing errors for each antenna at each timestep in radians. Has shape
(time, antenna, 2)
antenna_scaling :
numpy.ndarray
antenna scaling factor for each channel and each antenna. Has shape
(antenna, chan)
frequency :
numpy.ndarray
frequencies for each channel. Has shape
(chan,)
dtype :
numpy.dtype
, optionalNumpy dtype of result array. Should be float32 or float64. Defaults to float64
Returns: coords :
numpy.ndarray
coordinates of shape
(3, src, time, antenna, chan)
where each coordinate component represents l, m and frequency, respectively.- Rotating
-
africanus.rime.
beam_cube_dde
(beam, beam_lm_extents, beam_freq_map, lm, parallactic_angles, point_errors, antenna_scaling, frequency)[source]¶ Evaluates Direction Dependent Effects along a source’s path by interpolating the values of a complex beam cube at the source location.
Parameters: beam :
numpy.ndarray
Complex beam cube of shape
(beam_lw, beam_mh, beam_nud, corr, corr)
. beam_lw, beam_mh and beam_nud define the size of the cube in the l, m and frequency dimensions, respectively.beam_lm_extents :
numpy.ndarray
lm extents of the beam cube of shape
(2, 2)
.[[lower_l, upper_l], [lower_m, upper_m]]
.beam_freq_map :
numpy.ndarray
Beam frequency map of shape
(beam_nud,)
. This array is used to define interpolation along the(chan,)
dimension.lm :
numpy.ndarray
Source lm coordinates of shape
(source, 2)
. These coordinates are:- Scaled if the associated frequency lies outside the beam cube.
- Offset by pointing errors:
point_errors
- Rotated by parallactic angles:
parallactic_angles
. - Scaled by antenna scaling factors:
antenna_scaling
.
parallactic_angles :
numpy.ndarray
Parallactic angles of shape
(time, ant)
.point_errors :
numpy.ndarray
Pointing errors of shape
(time, ant, chan, 2)
.antenna_scaling :
numpy.ndarray
Antenna scaling factors of shape
(ant, chan, 2)
frequency :
numpy.ndarray
Frequencies of shape
(chan,)
.Returns: ddes :
numpy.ndarray
Direction Dependent Effects of shape
(source, time, ant, chan, corr, corr)
Notes
- Sources are clamped to the provided beam_lm_extents.
- Frequencies outside the cube (i.e. outside beam_freq_map) introduce linear scaling to the lm coordinates of a source.
-
africanus.rime.
zernike_dde
(coords, coeffs, noll_index)[source]¶ Computes Direction Dependent Effects by evaluating Zernicke Polynomials defined by coefficients
coeffs
and noll indexesnoll_index
at the specified coordinatescoords
.Decomposition of a voxel beam cube into Zernicke polynomial coefficients can be achieved through the use of the eidos package.
Parameters: coords :
numpy.ndarray
Float coordinates at which to evaluate the zernike polynomials. Has shape
(3, source, time, ant, chan)
. The three components in the first dimension represent l, m and frequency coordinates, respectively.coeffs :
numpy.ndarray
complex Zernicke polynomial coefficients. Has shape
(ant, chan, corr_1, ..., corr_n, poly)
wherepoly
is the number of polynomial coefficients andcorr_1, ..., corr_n
are a variable number of correlation dimensions.noll_index :
numpy.ndarray
Noll index associated with each polynomial coefficient. Has shape
(ant, chan, corr_1, ..., corr_n, poly)
.Returns: dde :
numpy.ndarray
complex values with shape
(source, time, ant, chan, corr_1, ..., corr_n)
-
africanus.rime.
wsclean_predict
(uvw, lm, source_type, flux, coeffs, log_poly, ref_freq, gauss_shape, frequency)[source]¶ Predict visibilities from a WSClean sky model.
Parameters: uvw :
numpy.ndarray
UVW coordinates of shape
(row, 3)
lm :
numpy.ndarray
Source LM coordinates of shape
(source, 2)
, in radians. Derived from theRa
andDec
fields.source_type :
numpy.ndarray
Strings defining the source type of shape
(source,)
. Should be either"POINT"
or"GAUSSIAN"
. Contains theType
field.flux :
numpy.ndarray
Source flux of shape
(source,)
. Contains theI
field.coeffs :
numpy.ndarray
Source Polynomial coefficients of shape
(source, coeffs)
. Contains theSpectralIndex
field.log_poly :
numpy.ndarray
Source polynomial type of shape
(source,)
. If True, logarithmic polynomials are used. If False, standard polynomials are used. Contains theLogarithmicSI
field.ref_freq :
numpy.ndarray
Source Reference frequency of shape
(source,)
. Contains theReferenceFrequency
field.gauss_shape :
numpy.ndarray
Gaussian shape parameters of shape
(source, 3)
used when the correspondingsource_type
is"GAUSSIAN"
. The 3 components should contain theMajorAxis
,MinorAxis
andOrientation
fields in radians, respectively.frequency :
numpy.ndarray
Frequency of shape
(chan,)
.Returns: visibilities :
numpy.ndarray
Complex visibilities of shape
(row, chan, 1)
Cuda¶
predict_vis (time_index, antenna1, antenna2) |
Multiply Jones terms together to form model visibilities according to the following formula: |
phase_delay (lm, uvw, frequency) |
Computes the phase delay (K) term: |
feed_rotation (parallactic_angles[, feed_type]) |
Computes the 2x2 feed rotation (L) matrix from the parallactic_angles . |
beam_cube_dde (beam, beam_lm_ext, …) |
Evaluates Direction Dependent Effects along a source’s path by interpolating the values of a complex beam cube at the source location. |
-
africanus.rime.cuda.
predict_vis
(time_index, antenna1, antenna2, dde1_jones=None, source_coh=None, dde2_jones=None, die1_jones=None, base_vis=None, die2_jones=None)[source]¶ Multiply Jones terms together to form model visibilities according to the following formula:
\[V_{pq} = G_{p} \left( B_{pq} + \sum_{s} E_{ps} X_{pqs} E_{qs}^H \right) G_{q}^H\]where for antenna \(p\) and \(q\), and source \(s\):
- \(B_{{pq}}\) represent base coherencies.
- \(E_{{ps}}\) represents Direction-Dependent Jones terms.
- \(X_{{pqs}}\) represents a coherency matrix (per-source).
- \(G_{{p}}\) represents Direction-Independent Jones terms.
Generally, \(E_{ps}\), \(G_{p}\), \(X_{pqs}\) should be formed by using the RIME API functions and combining them together with
einsum()
.Please read the Notes
Parameters: time_index :
cupy.ndarray
Time index used to look up the antenna Jones index for a particular baseline with shape
(row,)
. Obtainable viacp.unique(time, return_inverse=True)[1]
.antenna1 :
cupy.ndarray
Antenna 1 index used to look up the antenna Jones for a particular baseline. with shape
(row,)
.antenna2 :
cupy.ndarray
Antenna 2 index used to look up the antenna Jones for a particular baseline. with shape
(row,)
.dde1_jones :
cupy.ndarray
, optional\(E_{ps}\) Direction-Dependent Jones terms for the first antenna. shape
(source,time,ant,chan,corr_1,corr_2)
source_coh :
cupy.ndarray
, optional\(X_{pqs}\) Direction-Dependent Coherency matrix for the baseline. with shape
(source,row,chan,corr_1,corr_2)
dde2_jones :
cupy.ndarray
, optional\(E_{qs}\) Direction-Dependent Jones terms for the second antenna. This is usually the same array as
dde1_jones
as this preserves the symmetry of the RIME.predict_vis
will perform the conjugate transpose internally. shape(source,time,ant,chan,corr_1,corr_2)
die1_jones :
cupy.ndarray
, optional\(G_{ps}\) Direction-Independent Jones terms for the first antenna of the baseline. with shape
(time,ant,chan,corr_1,corr_2)
base_vis :
cupy.ndarray
, optional\(B_{pq}\) base coherencies, added to source coherency summation before multiplication with die1_jones and die2_jones. shape
(row,chan,corr_1,corr_2)
.die2_jones :
cupy.ndarray
, optional\(G_{ps}\) Direction-Independent Jones terms for the second antenna of the baseline. This is usually the same array as
die1_jones
as this preserves the symmetry of the RIME.predict_vis
will perform the conjugate transpose internally. shape(time,ant,chan,corr_1,corr_2)
Returns: visibilities :
cupy.ndarray
Model visibilities of shape
(row,chan,corr_1,corr_2)
Notes
- Direction-Dependent terms (dde{1,2}_jones) and Independent (die{1,2}_jones) are optional, but if one is present, the other must be present.
- The inputs to this function involve
row
,time
andant
(antenna) dimensions. - Each
row
is associated with a pair of antenna Jones matrices at a particular timestep via thetime_index
,antenna1
andantenna2
inputs. - The
row
dimension must be an increasing partial order in time.
-
africanus.rime.cuda.
phase_delay
(lm, uvw, frequency)[source]¶ Computes the phase delay (K) term:
\[ \begin{align}\begin{aligned}& {\Large e^{-2 \pi i (u l + v m + w (n - 1))} }\\& \textrm{where } n = \sqrt{1 - l^2 - m^2}\end{aligned}\end{align} \]Parameters: lm :
cupy.ndarray
LM coordinates of shape
(source, 2)
with L and M components in the last dimension.uvw :
cupy.ndarray
UVW coordinates of shape
(row, 3)
with U, V and W components in the last dimension.frequency :
cupy.ndarray
frequencies of shape
(chan,)
convention : {‘fourier’, ‘casa’}
Uses the \(e^{-2 \pi \mathit{i}}\) sign convention if
fourier
and \(e^{2 \pi \mathit{i}}\) ifcasa
.Returns: complex_phase :
cupy.ndarray
complex of shape
(source, row, chan)
Notes
Corresponds to the complex exponential of the Van Cittert-Zernike Theorem.
MeqTrees uses the CASA sign convention.
-
africanus.rime.cuda.
feed_rotation
(parallactic_angles, feed_type='linear')[source]¶ Computes the 2x2 feed rotation (L) matrix from the
parallactic_angles
.\[\begin{split}\textrm{linear} \begin{bmatrix} cos(pa) & sin(pa) \\ -sin(pa) & cos(pa) \end{bmatrix} \qquad \textrm{circular} \begin{bmatrix} e^{-i pa} & 0 \\ 0 & e^{i pa} \end{bmatrix}\end{split}\]Parameters: parallactic_angles :
cupy.ndarray
floating point parallactic angles. Of shape
(pa0, pa1, ..., pan)
.feed_type : {‘linear’, ‘circular’}
The type of feed
Returns: feed_matrix :
cupy.ndarray
Feed rotation matrix of shape
(pa0, pa1,...,pan,2,2)
-
africanus.rime.cuda.
beam_cube_dde
(beam, beam_lm_ext, beam_freq_map, lm, parangles, pointing_errors, antenna_scaling, frequencies)[source]¶ Evaluates Direction Dependent Effects along a source’s path by interpolating the values of a complex beam cube at the source location.
Parameters: beam :
cupy.ndarray
Complex beam cube of shape
(beam_lw, beam_mh, beam_nud, corr, corr)
. beam_lw, beam_mh and beam_nud define the size of the cube in the l, m and frequency dimensions, respectively.beam_lm_extents :
cupy.ndarray
lm extents of the beam cube of shape
(2, 2)
.[[lower_l, upper_l], [lower_m, upper_m]]
.beam_freq_map :
cupy.ndarray
Beam frequency map of shape
(beam_nud,)
. This array is used to define interpolation along the(chan,)
dimension.lm :
cupy.ndarray
Source lm coordinates of shape
(source, 2)
. These coordinates are:- Scaled if the associated frequency lies outside the beam cube.
- Offset by pointing errors:
point_errors
- Rotated by parallactic angles:
parallactic_angles
. - Scaled by antenna scaling factors:
antenna_scaling
.
parallactic_angles :
cupy.ndarray
Parallactic angles of shape
(time, ant)
.point_errors :
cupy.ndarray
Pointing errors of shape
(time, ant, chan, 2)
.antenna_scaling :
cupy.ndarray
Antenna scaling factors of shape
(ant, chan, 2)
frequency :
cupy.ndarray
Frequencies of shape
(chan,)
.Returns: ddes :
cupy.ndarray
Direction Dependent Effects of shape
(source, time, ant, chan, corr, corr)
Notes
- Sources are clamped to the provided beam_lm_extents.
- Frequencies outside the cube (i.e. outside beam_freq_map) introduce linear scaling to the lm coordinates of a source.
Dask¶
predict_vis (time_index, antenna1, antenna2) |
Multiply Jones terms together to form model visibilities according to the following formula: |
phase_delay (lm, uvw, frequency[, convention]) |
Computes the phase delay (K) term: |
parallactic_angles (times, antenna_positions, …) |
Computes parallactic angles per timestep for the given reference antenna position and field centre. |
feed_rotation (parallactic_angles, feed_type) |
Computes the 2x2 feed rotation (L) matrix from the parallactic_angles . |
transform_sources (lm, parallactic_angles, …) |
Creates beam sampling coordinates suitable for use in beam_cube_dde() by: |
beam_cube_dde (beam, beam_lm_extents, …) |
Evaluates Direction Dependent Effects along a source’s path by interpolating the values of a complex beam cube at the source location. |
zernike_dde (coords, coeffs, noll_index) |
Computes Direction Dependent Effects by evaluating Zernicke Polynomials defined by coefficients coeffs and noll indexes noll_index at the specified coordinates coords . |
wsclean_predict (uvw, lm, source_type, flux, …) |
Predict visibilities from a WSClean sky model. |
-
africanus.rime.dask.
predict_vis
(time_index, antenna1, antenna2, dde1_jones=None, source_coh=None, dde2_jones=None, die1_jones=None, base_vis=None, die2_jones=None, streams=None)[source]¶ Multiply Jones terms together to form model visibilities according to the following formula:
\[V_{pq} = G_{p} \left( B_{pq} + \sum_{s} E_{ps} X_{pqs} E_{qs}^H \right) G_{q}^H\]where for antenna \(p\) and \(q\), and source \(s\):
- \(B_{{pq}}\) represent base coherencies.
- \(E_{{ps}}\) represents Direction-Dependent Jones terms.
- \(X_{{pqs}}\) represents a coherency matrix (per-source).
- \(G_{{p}}\) represents Direction-Independent Jones terms.
Generally, \(E_{ps}\), \(G_{p}\), \(X_{pqs}\) should be formed by using the RIME API functions and combining them together with
einsum()
.Please read the Notes
Parameters: time_index :
dask.array.Array
Time index used to look up the antenna Jones index for a particular baseline with shape
(row,)
. Obtainable viatime.map_blocks(lambda a: np.unique(a, return_inverse=True)[1])
.antenna1 :
dask.array.Array
Antenna 1 index used to look up the antenna Jones for a particular baseline. with shape
(row,)
.antenna2 :
dask.array.Array
Antenna 2 index used to look up the antenna Jones for a particular baseline. with shape
(row,)
.dde1_jones :
dask.array.Array
, optional\(E_{ps}\) Direction-Dependent Jones terms for the first antenna. shape
(source,time,ant,chan,corr_1,corr_2)
source_coh :
dask.array.Array
, optional\(X_{pqs}\) Direction-Dependent Coherency matrix for the baseline. with shape
(source,row,chan,corr_1,corr_2)
dde2_jones :
dask.array.Array
, optional\(E_{qs}\) Direction-Dependent Jones terms for the second antenna. This is usually the same array as
dde1_jones
as this preserves the symmetry of the RIME.predict_vis
will perform the conjugate transpose internally. shape(source,time,ant,chan,corr_1,corr_2)
die1_jones :
dask.array.Array
, optional\(G_{ps}\) Direction-Independent Jones terms for the first antenna of the baseline. with shape
(time,ant,chan,corr_1,corr_2)
base_vis :
dask.array.Array
, optional\(B_{pq}\) base coherencies, added to source coherency summation before multiplication with die1_jones and die2_jones. shape
(row,chan,corr_1,corr_2)
.die2_jones :
dask.array.Array
, optional\(G_{ps}\) Direction-Independent Jones terms for the second antenna of the baseline. This is usually the same array as
die1_jones
as this preserves the symmetry of the RIME.predict_vis
will perform the conjugate transpose internally. shape(time,ant,chan,corr_1,corr_2)
streams : {False, True}
If
True
the coherencies are serially summed in a linear chain. IfFalse
, dask uses a tree style reduction algorithm.Returns: visibilities :
dask.array.Array
Model visibilities of shape
(row,chan,corr_1,corr_2)
Notes
Direction-Dependent terms (dde{1,2}_jones) and Independent (die{1,2}_jones) are optional, but if one is present, the other must be present.
The inputs to this function involve
row
,time
andant
(antenna) dimensions.Each
row
is associated with a pair of antenna Jones matrices at a particular timestep via thetime_index
,antenna1
andantenna2
inputs.The
row
dimension must be an increasing partial order in time.The
ant
dimension should only contain a single chunk equal to the number of antenna. Since eachrow
can contain any antenna, random access must be preserved along this dimension.The chunks in the
row
andtime
dimension must align. This subtle point must be understood otherwise invalid results will be produced by the chunking scheme. In the example below we have four unique time indices[0,1,2,3]
, and four unique antenna[0,1,2,3]
indexing10
rows.# Row indices into the time/antenna indexed arrays time_idx = np.asarray([0,0,1,1,2,2,2,2,3,3]) ant1 = np.asarray( [0,0,0,0,1,1,1,2,2,3] ant2 = np.asarray( [0,1,2,3,1,2,3,2,3,3])
A reasonable chunking scheme for the
row
andtime
dimension would be(4,4,2)
and(2,1,1)
respectively. Another way of explaining this is that the first four rows contain two unique timesteps, the second four rows contain one unique timestep and the last two rows contain one unique timestep.Some rules of thumb:
The number chunks in
row
andtime
must match although the individual chunk sizes need not.Unique timesteps should not be split across row chunks.
For a Measurement Set whose rows are ordered on the
TIME
column, the following is a good way of obtaining the row chunking strategy:import numpy as np import pyrap.tables as pt ms = pt.table("data.ms") times = ms.getcol("TIME") unique_times, chunks = np.unique(times, return_counts=True)
Use
aggregate_chunks()
to aggregate multiplerow
andtime
chunks into chunks large enough such that functions operating on the resulting data can drop the GIL and spend time processing the data. Expanding the previous example:# Aggregate row utimes = unique_times.size # Single chunk for each unique time time_chunks = (1,)*utimes # Aggregate row chunks into chunks <= 10000 aggregate_chunks((chunks, time_chunks), (10000, utimes))
-
africanus.rime.dask.
phase_delay
(lm, uvw, frequency, convention='fourier')[source]¶ Computes the phase delay (K) term:
\[ \begin{align}\begin{aligned}& {\Large e^{-2 \pi i (u l + v m + w (n - 1))} }\\& \textrm{where } n = \sqrt{1 - l^2 - m^2}\end{aligned}\end{align} \]Parameters: lm :
dask.array.Array
LM coordinates of shape
(source, 2)
with L and M components in the last dimension.uvw :
dask.array.Array
UVW coordinates of shape
(row, 3)
with U, V and W components in the last dimension.frequency :
dask.array.Array
frequencies of shape
(chan,)
convention : {‘fourier’, ‘casa’}
Uses the \(e^{-2 \pi \mathit{i}}\) sign convention if
fourier
and \(e^{2 \pi \mathit{i}}\) ifcasa
.Returns: complex_phase :
dask.array.Array
complex of shape
(source, row, chan)
Notes
Corresponds to the complex exponential of the Van Cittert-Zernike Theorem.
MeqTrees uses the CASA sign convention.
-
africanus.rime.dask.
parallactic_angles
(times, antenna_positions, field_centre, **kwargs)[source]¶ Computes parallactic angles per timestep for the given reference antenna position and field centre.
Parameters: times :
dask.array.Array
Array of Mean Julian Date times in seconds with shape
(time,)
,antenna_positions :
dask.array.Array
Antenna positions of shape
(ant, 3)
in metres in the ITRF frame.field_centre :
dask.array.Array
Field centre of shape
(2,)
in radiansbackend : {‘casa’, ‘test’}, optional
Backend to use for calculating the parallactic angles.
casa
defers to an implementation depending onpython-casacore
. This backend should be used by default.test
creates parallactic angles by multiplying thetimes
andantenna_position
arrays. It exist solely for testing.
Returns: parallactic_angles :
dask.array.Array
Parallactic angles of shape
(time,ant)
-
africanus.rime.dask.
feed_rotation
(parallactic_angles, feed_type)[source]¶ Computes the 2x2 feed rotation (L) matrix from the
parallactic_angles
.\[\begin{split}\textrm{linear} \begin{bmatrix} cos(pa) & sin(pa) \\ -sin(pa) & cos(pa) \end{bmatrix} \qquad \textrm{circular} \begin{bmatrix} e^{-i pa} & 0 \\ 0 & e^{i pa} \end{bmatrix}\end{split}\]Parameters: parallactic_angles :
numpy.ndarray
floating point parallactic angles. Of shape
(pa0, pa1, ..., pan)
.feed_type : {‘linear’, ‘circular’}
The type of feed
Returns: feed_matrix :
numpy.ndarray
Feed rotation matrix of shape
(pa0, pa1,...,pan,2,2)
-
africanus.rime.dask.
transform_sources
(lm, parallactic_angles, pointing_errors, antenna_scaling, frequency, dtype=None)[source]¶ Creates beam sampling coordinates suitable for use in
beam_cube_dde()
by:- Rotating
lm
coordinates by theparallactic_angles
- Adding
pointing_errors
- Scaling by
antenna_scaling
Parameters: lm :
dask.array.Array
LM coordinates of shape
(src,2)
in radians offset from the phase centre.parallactic_angles :
dask.array.Array
parallactic angles of shape
(time, antenna)
in radians.pointing_errors :
dask.array.Array
LM pointing errors for each antenna at each timestep in radians. Has shape
(time, antenna, 2)
antenna_scaling :
dask.array.Array
antenna scaling factor for each channel and each antenna. Has shape
(antenna, chan)
frequency :
dask.array.Array
frequencies for each channel. Has shape
(chan,)
dtype :
numpy.dtype
, optionalNumpy dtype of result array. Should be float32 or float64. Defaults to float64
Returns: coords :
dask.array.Array
coordinates of shape
(3, src, time, antenna, chan)
where each coordinate component represents l, m and frequency, respectively.- Rotating
-
africanus.rime.dask.
beam_cube_dde
(beam, beam_lm_extents, beam_freq_map, lm, parallactic_angles, point_errors, antenna_scaling, frequencies)[source]¶ Evaluates Direction Dependent Effects along a source’s path by interpolating the values of a complex beam cube at the source location.
Parameters: beam :
dask.array.Array
Complex beam cube of shape
(beam_lw, beam_mh, beam_nud, corr, corr)
. beam_lw, beam_mh and beam_nud define the size of the cube in the l, m and frequency dimensions, respectively.beam_lm_extents :
dask.array.Array
lm extents of the beam cube of shape
(2, 2)
.[[lower_l, upper_l], [lower_m, upper_m]]
.beam_freq_map :
dask.array.Array
Beam frequency map of shape
(beam_nud,)
. This array is used to define interpolation along the(chan,)
dimension.lm :
dask.array.Array
Source lm coordinates of shape
(source, 2)
. These coordinates are:- Scaled if the associated frequency lies outside the beam cube.
- Offset by pointing errors:
point_errors
- Rotated by parallactic angles:
parallactic_angles
. - Scaled by antenna scaling factors:
antenna_scaling
.
parallactic_angles :
dask.array.Array
Parallactic angles of shape
(time, ant)
.point_errors :
dask.array.Array
Pointing errors of shape
(time, ant, chan, 2)
.antenna_scaling :
dask.array.Array
Antenna scaling factors of shape
(ant, chan, 2)
frequency :
dask.array.Array
Frequencies of shape
(chan,)
.Returns: ddes :
dask.array.Array
Direction Dependent Effects of shape
(source, time, ant, chan, corr, corr)
Notes
- Sources are clamped to the provided beam_lm_extents.
- Frequencies outside the cube (i.e. outside beam_freq_map) introduce linear scaling to the lm coordinates of a source.
-
africanus.rime.dask.
zernike_dde
(coords, coeffs, noll_index)[source]¶ Computes Direction Dependent Effects by evaluating Zernicke Polynomials defined by coefficients
coeffs
and noll indexesnoll_index
at the specified coordinatescoords
.Decomposition of a voxel beam cube into Zernicke polynomial coefficients can be achieved through the use of the eidos package.
Parameters: coords :
dask.array.Array
Float coordinates at which to evaluate the zernike polynomials. Has shape
(3, source, time, ant, chan)
. The three components in the first dimension represent l, m and frequency coordinates, respectively.coeffs :
dask.array.Array
complex Zernicke polynomial coefficients. Has shape
(ant, chan, corr_1, ..., corr_n, poly)
wherepoly
is the number of polynomial coefficients andcorr_1, ..., corr_n
are a variable number of correlation dimensions.noll_index :
dask.array.Array
Noll index associated with each polynomial coefficient. Has shape
(ant, chan, corr_1, ..., corr_n, poly)
.Returns: dde :
dask.array.Array
complex values with shape
(source, time, ant, chan, corr_1, ..., corr_n)
-
africanus.rime.dask.
wsclean_predict
(uvw, lm, source_type, flux, coeffs, log_poly, ref_freq, gauss_shape, frequency)[source]¶ Predict visibilities from a WSClean sky model.
Parameters: uvw :
dask.array.Array
UVW coordinates of shape
(row, 3)
lm :
dask.array.Array
Source LM coordinates of shape
(source, 2)
, in radians. Derived from theRa
andDec
fields.source_type :
dask.array.Array
Strings defining the source type of shape
(source,)
. Should be either"POINT"
or"GAUSSIAN"
. Contains theType
field.flux :
dask.array.Array
Source flux of shape
(source,)
. Contains theI
field.coeffs :
dask.array.Array
Source Polynomial coefficients of shape
(source, coeffs)
. Contains theSpectralIndex
field.log_poly :
dask.array.Array
Source polynomial type of shape
(source,)
. If True, logarithmic polynomials are used. If False, standard polynomials are used. Contains theLogarithmicSI
field.ref_freq :
dask.array.Array
Source Reference frequency of shape
(source,)
. Contains theReferenceFrequency
field.gauss_shape :
dask.array.Array
Gaussian shape parameters of shape
(source, 3)
used when the correspondingsource_type
is"GAUSSIAN"
. The 3 components should contain theMajorAxis
,MinorAxis
andOrientation
fields in radians, respectively.frequency :
dask.array.Array
Frequency of shape
(chan,)
.Returns: visibilities :
dask.array.Array
Complex visibilities of shape
(row, chan, 1)