#include <stdio.h>
#include <stdlib.h>
#include <unistd.h>
#include <time.h>
#include <math.h>
#include <string.h>
#include <mwalib.h>
#include "vcsbeam.h"
#include "gpu_macros.h"

Include dependency graph for jones.c:

Functions
void	vmSetPolIdxLists (vcsbeam_context *vm)
	Creates arrays of indexes for antennas and polarisation according to the ordering used in legacy VCS observations. More...

void	vmCalcJ (vcsbeam_context *vm)
	Computes the (inverse) Jones matrix \({\bf J}^{-1} = \left({\bf D}{\bf B}\right)^{-1}\). More...

void	vmCalcJonesAndDelays (vcsbeam_context vm, double ras_hours, double decs_degs, beam_geom beam_geom_vals)
	Wrapper function for vmCalcPhi(), vmCalcB(), and vmCalcJ(). More...

void	cp2x2 (gpuDoubleComplex Min, gpuDoubleComplex Mout)
	Copies a \(2\times2\) complex-valued matrix. More...

gpuDoubleComplex	reciprocal_complex (gpuDoubleComplex z)
	Computes the reciprocal of a complex number. More...

gpuDoubleComplex	negate_complex (gpuDoubleComplex z)
	Computes the negative of a complex number. More...

void	inv2x2 (gpuDoubleComplex Min, gpuDoubleComplex Mout)
	Computes the inverse of a \(2\times2\) complex-valued matrix. More...

void	inv2x2d (double Min, double Mout)
	Computes the inverse of a \(2\times2\) real-valued matrix. More...

void	inv2x2S (gpuDoubleComplex Min, gpuDoubleComplex Mout)
	Computes the inverse of a \(2\times2\) complex-valued matrix. More...

void	mult2x2d (gpuDoubleComplex M1, gpuDoubleComplex M2, gpuDoubleComplex *Mout)
	Performs a matrix multiplication of two \(2\times2\) complex-valued matrices. More...

void	mult2x2d_RxC (double M1, gpuDoubleComplex M2, gpuDoubleComplex *Mout)
	Performs a matrix multiplication of a \(2\times2\) real-valued matrix (on the left) with a \(2\times2\) complex-valued matrix (on the right). More...

void	mult2x2d_CxR (gpuDoubleComplex M1, double M2, gpuDoubleComplex *Mout)
	Performs a matrix multiplication of a \(2\times2\) complex-valued matrix (on the left) with a \(2\times2\) real-valued matrix (on the right). More...

void	conj2x2 (gpuDoubleComplex M, gpuDoubleComplex Mout)
	Computes the conjugate of a \(2\times2\) complex-valued matrix. More...

double	norm2x2 (gpuDoubleComplex M, gpuDoubleComplex Mout)
	Normalises a \(2\times2\) matrix via the Frobenius norm. More...

void	reverse2x2 (gpuDoubleComplex M, gpuDoubleComplex Mout)
	Reverses the elements of a \(2\times2\) matrix. More...

void	swaprows2x2 (gpuDoubleComplex M, gpuDoubleComplex Mout)
	Swaps the rows of a \(2\times2\) matrix. More...

void	swapcols2x2 (gpuDoubleComplex M, gpuDoubleComplex Mout)
	Swaps the columns of a \(2\times2\) matrix. More...

bool	is2x2zero (gpuDoubleComplex *M)
	Returns true if and only if all elements of a \(2\times2\) complex-valued matrix are identically zero. More...

void	calc_hermitian (gpuDoubleComplex M, gpuDoubleComplex H)
	Computes the Hermitian (i.e. More...

void	calc_coherency_matrix (gpuDoubleComplex M, gpuDoubleComplex C)
	Computes the coherency matrix of a \(2\times2\) complex-valued matrix. More...

void	fprintf_complex_matrix (FILE fout, gpuDoubleComplex M)
	Prints a \(2\times2\) complex-valued matrix. More...

Function Documentation

◆ calc_coherency_matrix()

void calc_coherency_matrix	(	gpuDoubleComplex *	M,
		gpuDoubleComplex *	C
	)

Computes the coherency matrix of a \(2\times2\) complex-valued matrix.

Parameters

	M	The input matrix, \({\bf M}\).
[out]	C	The output matrix, \({\bf C} = {\bf M}{\bf M}^\dagger\).

It is safe for M and C to point to the same matrix.

◆ calc_hermitian()

void calc_hermitian	(	gpuDoubleComplex *	M,
		gpuDoubleComplex *	H
	)

Computes the Hermitian (i.e.

conjugate transpose) of a \(2\times2\) complex-valued matrix.

Parameters

	M	The input matrix, \({\bf M}\).
[out]	H	The output matrix, \({\bf H} = {\bf M}^\dagger\).

It is safe for M and H to point to the same matrix.

◆ conj2x2()

void conj2x2	(	gpuDoubleComplex *	M,
		gpuDoubleComplex *	Mout
	)

Computes the conjugate of a \(2\times2\) complex-valued matrix.

Parameters

	M	The input matrix, \({\bf M}\).
[out]	Mout	The output matrix, \({\bf M}_\text{out} = {\bf M}^\ast\).

It is safe for M and Mout to point to the same matrix.

◆ cp2x2()

void cp2x2	(	gpuDoubleComplex *	Min,
		gpuDoubleComplex *	Mout
	)

Copies a \(2\times2\) complex-valued matrix.

Parameters

	Min	The source matrix, \({\bf M}_\text{in}\).
[out]	Mout	The destination matrix, \({\bf M}_\text{out} = {\bf M}_\text{in}\).

◆ fprintf_complex_matrix()

void fprintf_complex_matrix	(	FILE *	fout,
		gpuDoubleComplex *	M
	)

Prints a \(2\times2\) complex-valued matrix.

Parameters

fout	The file stream to output to
M	The matrix to print

The output format is

[ a+b*i, c+d*i; e+f*i, g+h*i ]

and a newline is printed at the end.

◆ inv2x2()

void inv2x2	(	gpuDoubleComplex *	Min,
		gpuDoubleComplex *	Mout
	)

Computes the inverse of a \(2\times2\) complex-valued matrix.

Parameters

	Min	The input matrix, \({\bf M}_\text{in}\).
[out]	Mout	The output matrix, \({\bf M}_\text{out} = {\bf M}_\text{in}^{-1}\)

◆ inv2x2d()

void inv2x2d	(	double *	Min,
		double *	Mout
	)

Computes the inverse of a \(2\times2\) real-valued matrix.

Parameters

	Min	The input matrix, \({\bf M}_\text{in}\)
[out]	Mout	The output matrix, \({\bf M}_\text{out} = {\bf M}_\text{in}^{-1}\).

◆ inv2x2S()

void inv2x2S	(	gpuDoubleComplex *	Min,
		gpuDoubleComplex *	Mout
	)

Computes the inverse of a \(2\times2\) complex-valued matrix.

Parameters

	Min	The input matrix, \({\bf M}_\text{in}\)
[out]	Mout	The output matrix, \({\bf M}_\text{out} = {\bf M}_\text{in}^{-1}\)

See also: inv2x2()

Todo:: Does this do the same thing as inv2x2()? Check the code, and delete if not needed.

◆ is2x2zero()

bool is2x2zero ( gpuDoubleComplex * M )

Returns true if and only if all elements of a \(2\times2\) complex-valued matrix are identically zero.

Parameters

M	The matrix to be tested, \({\bf M}\).

◆ mult2x2d()

void mult2x2d	(	gpuDoubleComplex *	M1,
		gpuDoubleComplex *	M2,
		gpuDoubleComplex *	Mout
	)

Performs a matrix multiplication of two \(2\times2\) complex-valued matrices.

Parameters

	M1	The first input matrix, \({\bf M}_1\)
	M2	The second input matrix, \({\bf M}_2\)
[out]	Mout	The output matrix, \({\bf M}_\text{out} = {\bf M}_1 \times {\bf M}_2\).

◆ mult2x2d_CxR()

void mult2x2d_CxR	(	gpuDoubleComplex *	M1,
		double *	M2,
		gpuDoubleComplex *	Mout
	)

Performs a matrix multiplication of a \(2\times2\) complex-valued matrix (on the left) with a \(2\times2\) real-valued matrix (on the right).

Parameters

	M1	The first input matrix (complex-valued), \({\bf M}_1\)
	M2	The second input matrix (real-valued), \({\bf M}_2\)
[out]	Mout	The output matrix (complex-valued), \({\bf M}_\text{out} = {\bf M}_1 \times {\bf M}_2\).

◆ mult2x2d_RxC()

void mult2x2d_RxC	(	double *	M1,
		gpuDoubleComplex *	M2,
		gpuDoubleComplex *	Mout
	)

Performs a matrix multiplication of a \(2\times2\) real-valued matrix (on the left) with a \(2\times2\) complex-valued matrix (on the right).

Parameters

	M1	The first input matrix (real-valued), \({\bf M}_1\)
	M2	The second input matrix (complex-valued), \({\bf M}_2\)
[out]	Mout	The output matrix (complex-valued), \({\bf M}_\text{out} = {\bf M}_1 \times {\bf M}_2\).

◆ negate_complex()

gpuDoubleComplex negate_complex ( gpuDoubleComplex z )

Computes the negative of a complex number.

Parameters

z	A complex number

Returns: \(-z\)

◆ norm2x2()

double norm2x2	(	gpuDoubleComplex *	M,
		gpuDoubleComplex *	Mout
	)

Normalises a \(2\times2\) matrix via the Frobenius norm.

Parameters

	M	The input matrix, \({\bf M}\).
[out]	Mout	The output matrix, \({\bf M}_\text{out} = \|{\bf M}_\text{in}\|\).

It is safe for M and Mout to point to the same matrix.

◆ reciprocal_complex()

gpuDoubleComplex reciprocal_complex ( gpuDoubleComplex z )

Computes the reciprocal of a complex number.

Parameters

z	A complex number

Returns: \(\dfrac{1}{z}\)

◆ reverse2x2()

void reverse2x2	(	gpuDoubleComplex *	M,
		gpuDoubleComplex *	Mout
	)

Reverses the elements of a \(2\times2\) matrix.

Parameters

	M	The input matrix, \({\bf M}\).
[out]	Mout	The output matrix, \({\bf M}_\text{out}\).

\[\begin{bmatrix} 0 & 1 \\ 2 & 3\end{bmatrix} \rightarrow \begin{bmatrix} 3 & 2 \\ 1 & 0 \end{bmatrix}\]

It is safe for M and Mout to point to the same matrix.

◆ swapcols2x2()

void swapcols2x2	(	gpuDoubleComplex *	M,
		gpuDoubleComplex *	Mout
	)

Swaps the columns of a \(2\times2\) matrix.

Parameters

	M	The input matrix, \({\bf M}\).
[out]	Mout	The output matrix, \({\bf M}_\text{out}\).

\[\begin{bmatrix} 0 & 1 \\ 2 & 3\end{bmatrix} \rightarrow \begin{bmatrix} 1 & 0 \\ 3 & 2 \end{bmatrix}\]

It is safe for M and Mout to point to the same matrix.

◆ swaprows2x2()

void swaprows2x2	(	gpuDoubleComplex *	M,
		gpuDoubleComplex *	Mout
	)

Swaps the rows of a \(2\times2\) matrix.

Parameters

	M	The input matrix, \({\bf M}\).
[out]	Mout	The output matrix, \({\bf M}_\text{out}\).

\[\begin{bmatrix} 0 & 1 \\ 2 & 3\end{bmatrix} \rightarrow \begin{bmatrix} 2 & 3 \\ 0 & 1 \end{bmatrix}\]

It is safe for M and Mout to point to the same matrix.

◆ vmCalcJ()

void vmCalcJ ( vcsbeam_context * vm )

Computes the (inverse) Jones matrix \({\bf J}^{-1} = \left({\bf D}{\bf B}\right)^{-1}\).

\({\bf J}^{-1}\) is the Jones matrix that converts from instrumental voltages to sky voltages, according to

\[ \begin{aligned} {\bf e} &= {\bf J}^{-1} {\bf v}, \\ \begin{bmatrix} e_x \\ e_y \end{bmatrix} &= \begin{bmatrix} J_{xq} & J_{xp} \\ J_{yq} & J_{yp} \end{bmatrix} \begin{bmatrix} v_q \\ v_p \end{bmatrix}. \end{aligned} \]

It is computed by multiplying the instrumental calibration Jones matrix with the beam model matrix, and taking the inverse

\[{\bf J}^{-1} = \left({\bf D}{\bf B}\right)^{-1},\]

where \({\bf D}\) is taken from vm→D, and \({\bf B}\) from vm→pb.B. The result is stored in vm→J.

\({\bf J}^{-1}\) has dimensions \((N_a \times N_f \times N_p \times N_p)\), where

\(N_a\) is the number of antennas (vm→obs_metadata→num_ants)
\(N_f\) is the number of frequencies (vm→vm→nfine_chan)
\(N_p\) is the number of polarisation (vm→obs_metadata→num_ant_pols)

The antenna ordering is given by the "Antenna" keyword in the observation's metafits file. The frequencies are ordered from lowest to highest.

Todo:: Issue #9

◆ vmCalcJonesAndDelays()

void vmCalcJonesAndDelays	(	vcsbeam_context *	vm,
		double *	ras_hours,
		double *	decs_degs,
		beam_geom *	beam_geom_vals
	)

Wrapper function for vmCalcPhi(), vmCalcB(), and vmCalcJ().

Parameters

vm	The VCSBeam context struct
ras_hours	An array of right ascensions, in decimal hours, one for each pointing.
decs_degs	An array of declinations, in decimal degrees, one for each pointing.
beam_geom_vals	An array of beam geometries, one for each pointing.

For each pointing, the quantities \(e^{i\varphi}\), \({\bf B}\), and \({\bf J}^{-1}\) are calculated. ( \({\bf D}\) is not recalculated, but is used in the calculation of \({\bf J}^{-1}\)).

◆ vmSetPolIdxLists()

void vmSetPolIdxLists ( vcsbeam_context * vm )

Creates arrays of indexes for antennas and polarisation according to the ordering used in legacy VCS observations.

Parameters

vm	The VCSBeam context struct

The index arrays are needed for pairing the antenna metadata with the voltage samples, which historically were ordered in a different way. Fortunately, the metafits file contains the VCS ordering in the "VCSOrder" field. This function converts these indexes into a lookup array that can be conveniently passed to GPU kernels at runtime. A separate array is made for each polarisation, with the indexes into the arrays being the metafits "Antenna" number, and the value stored in that array position being the "VCSOrder".

The metafits information is drawn from vm→obs_metadata, and the indexes are saved in vm→polP_idxs and vm→polQ_idxs.

No memory is allocated by this function. It is assumed that all arrays have been previously allocated, and are sufficiently large.

Functions

Function Documentation

◆ calc_coherency_matrix()

◆ calc_hermitian()

◆ conj2x2()

◆ cp2x2()

◆ fprintf_complex_matrix()

◆ inv2x2()

◆ inv2x2d()

◆ inv2x2S()

◆ is2x2zero()

◆ mult2x2d()

◆ mult2x2d_CxR()

◆ mult2x2d_RxC()

◆ negate_complex()

◆ norm2x2()

◆ reciprocal_complex()

◆ reverse2x2()

◆ swapcols2x2()

◆ swaprows2x2()

◆ vmCalcJ()

◆ vmCalcJonesAndDelays()

◆ vmSetPolIdxLists()