6. Python API reference

The documentation of each function that is reported in this page can also be directly accessed from Python with, e.g.:

from galario.double import sampleImage

help(sampleImage)

6.1. Computing synthetic visibilities

To compute the synthetic visibilities of a model use the sampleImage() and sampleProfile() functions.

galario.double.sampleImage()

Compute the synthetic visibilities of a model image at the specified (u, v) locations.

The 2D surface brightness in image is Fourier transformed and sampled in the (u, v) locations given in the u and v arrays.

Typical call signature:

vis = sampleImage(image, dxy, u, v, dRA=0, dDec=0, PA=0, check=False, origin='upper')

Parameters:

• image (2D array_like, float) – Square matrix of shape (nxy, nxy) containing the 2D surface brightness of the model. Assume the x-axis (R.A.) increases from right (West) to left (East) and the y-axis (Dec.) increases from bottom (South) to top (North). nxy must be even. units: Jy/pixel

• dxy (float) – Size of the image cell, assumed equal in both x and y direction. units: rad

• u (array_like, float) – u coordinate of the visibility points where the FT has to be sampled. units: wavelength

• v (array_like, float) – v coordinate of the visibility points where the FT has to be sampled. The length of v must be equal to the length of u. units: wavelength

• dRA (float, optional) – R.A. offset w.r.t. the phase center by which the image is translated. If dRA > 0 translate the image towards the left (East). Default is 0. units: rad

• dDec (float, optional) – Dec. offset w.r.t. the phase center by which the image is translated. If dDec > 0 translate the image towards the top (North). Default is 0. units: rad

• PA (float, optional) – Position Angle, defined East of North. Default is 0. units: rad

• check (bool, optional) – If True, check whether image and dxy satisfy Nyquist criterion for computing the synthetic visibilities in the (u, v) locations provided. Additionally check that the (u, v) points fall in the image to avoid segmentation violations. Default is False since the check might take time. For executions where speed is important, set to False.

• origin (['upper' | 'lower'], optional) – Set the [0,0] pixel index of the matrix in the upper left or lower left corner of the axes. It follows the same convention as in matplotlib matshow and imshow commands. Declination axis and the matrix y axis are parallel for origin='lower', anti-parallel for origin='upper'. The central pixel corresponding to the (RA, Dec) = (0, 0) is always [Nxy/2, Nxy/2]. For more details see the Technical Requirements page in the online docs.

Returns:

vis – Synthetic visibilities sampled in the (u, v) locations given in u and v. units: Jy

Return type:

array_like, complex

galario.double.sampleProfile()

Compute the synthetic visibilities of a model with an axisymmetric brightness profile.

The brightness profile intensity is used to build a 2D image of the model, which is then Fourier transformed and sampled in the (u, v) locations given in the u and v arrays.

The image is created as in sweep() assuming that the x-axis (R.A.) increases from right (West) to left (East) and the y-axis (Dec.) increases from bottom (South) to top (North).

Typical call signature:

vis = sampleProfile(intensity, Rmin, dR, nxy, dxy, u, v,
                    dRA=0, dDec=0, PA=0, inc=0, check=False)

Parameters:

• intensity ((M,) array_like, float) – Array containing the radial brightness profile of the model. The profile is assumed to be sampled on a linear radial grid starting at Rmin with spacing dR. units: Jy/sr

• Rmin (float) – Inner edge of the radial grid, i.e. the radius where the brightness is intensity[0]. units: rad

• dR (float) – Size of the cell of the radial grid, assumed linear. units: rad

• nxy (int) – Side of the square model image, which is internally computed. units: pixel

• dxy (float) – Size of the image cell, assumed equal in both x and y direction. units: rad

• u (array_like, float) – u coordinate of the visibility points where the FT has to be sampled. units: wavelength

• v (array_like, float) – v coordinate of the visibility points where the FT has to be sampled. The length of v must be equal to the length of u. units: wavelength

• dRA (float, optional) – R.A. offset w.r.t. the phase center by which the image is translated. If dRA > 0 translate the image towards the left (East). Default is 0. units: rad

• dDec (float, optional) – Dec. offset w.r.t. the phase center by which the image is translated. If dDec > 0 translate the image towards the top (North). Default is 0. units: rad

• PA (float, optional) – Position Angle, defined East of North. Default is 0. units: rad

• inc (float, optional) – Inclination of the image plane along a North-South (top-bottom) axis. If inc=0. the image is face-on; if inc=90. the image is edge-on. units: rad

• check (bool, optional) – If True, check whether image and dxy satisfy Nyquist criterion for computing the synthetic visibilities in the (u, v) locations provided. Additionally check that the (u, v) points fall in the image to avoid segmentation violations. Default is False since the check might take time. For executions where speed is important, set to False.

Returns:

vis – Synthetic visibilities sampled in the (u, v) locations given in u and v. units: Jy

Return type:

array_like, complex

See also

sweep()

Note

The translation by (dRA, dDec) and the rotation by PA of the model image are optimized for speed: they are not obtained with extra interpolations of the model image, but rather applied in the Fourier plane.

The offset is achieved by applying a complex phase to the sampled visibilities. To rotation is achieved by internally rotating the (u, v) locations by -PA.

Suggestion: We recommend starting with check set to True to ensure the results obtained are correct. Once a combination of matrix size and dxy for the given data has been found, uvcheck can be safely set to False.

6.2. Computing chi squared directly

To compute the \(\chi^2\) of a model use the chi2Image() and chi2Profile() functions.

galario.double.chi2Image()

Compute the chi square of a model image given the observed visibilities.

The chi square is computed from the observed and synthetic visibilities as:

\[\chi^2 = \sum_{j=1}^N w_j * [(Re V_{obs\ j}-Re V_{mod\ j})^2 + (Im V_{obs\ j}-Im V_{mod\ j})^2]\]

where \(V_{mod}\) are the synthetic visibilities, which are computed internally as in sampleImage().

Typical call signature:

chi2 = chi2Image(image, dxy, u, v, vis_obs_re, vis_obs_im, vis_obs_w,
                 dRA=0, dDec=0, PA=0, check=False, origin='upper')

Parameters:

• image (2D array_like, float) – Square matrix of shape (nxy, nxy) containing the 2D surface brightness of the model. Assume the x-axis (R.A.) increases from right (West) to left (East) and the y-axis (Dec.) increases from bottom (South) to top (North). nxy must be even. units: Jy/pixel

• dxy (float) – Size of the image cell, assumed equal in both x and y direction. units: rad

• u (array_like, float) – u coordinate of the visibility points where the FT has to be sampled. units: wavelength

• v (array_like, float) – v coordinate of the visibility points where the FT has to be sampled. The length of v must be equal to the length of u. units: wavelength

• vis_obs_re (array_like, float) – Real part of the observed visibilities. units: Jy

• vis_obs_im (array_like, float) – Imaginary part of the observed visibilities. The length of vis_obs_im must be equal to the length of vis_obs_re. units: Jy

• vis_obs_w (array_like, float) – Weight associated to the observed visibilities. The length of vis_obs_w must be equal to the length of vis_obs_re. units:

• dRA (float, optional) – R.A. offset w.r.t. the phase center by which the image is translated. If dRA > 0 translate the image towards the left (East). Default is 0. units: rad

• dDec (float, optional) – Dec. offset w.r.t. the phase center by which the image is translated. If dDec > 0 translate the image towards the top (North). Default is 0. units: rad

• PA (float, optional) – Position Angle, defined East of North. Default is 0. units: rad

• check (bool, optional) – If True, check whether image and dxy satisfy Nyquist criterion for computing the synthetic visibilities in the (u, v) locations provided. Additionally check that the (u, v) points fall in the image to avoid segmentation violations. Default is False since the check might take time. For executions where speed is important, set to False.

• origin (['upper' | 'lower'], optional) – Set the [0,0] pixel index of the matrix in the upper left or lower left corner of the axes. It follows the same convention as in matplotlib matshow and imshow commands. Declination axis and the matrix y axis are parallel for origin='lower', anti-parallel for origin='upper'. The central pixel corresponding to the (RA, Dec) = (0, 0) is always [Nxy/2, Nxy/2]. For more details see the Technical Requirements page in the online docs.

Returns:

chi2 – The chi square, not normalized.

Return type:

float

See also

sampleImage()

galario.double.chi2Profile()

Compute the chi square of a model with an axisymmetric brightness profile given the observed visibilities.

The image is created from the intensity profile as in sweep(). The chi square is computed from the observed and synthetic visibilities as:

\[\chi^2 = \sum_{j=1}^N w_j * [(Re V_{obs\ j}-Re V_{mod\ j})^2 + (Im V_{obs\ j}-Im V_{mod\ j})^2]\]

where \(V_{mod}\) are the synthetic visibilities, which are computed internally as in sampleProfile().

Typical call signature:

chi2 = chi2Profile(intensity, Rmin, dR, nxy, dxy, u, v, vis_obs_re, vis_obs_im, vis_obs_w,
                   dRA=0, dDec=0, PA=0, inc=0, check=False)

Parameters:

• intensity (array_like, float) – Array containing the radial brightness profile of the model. The profile is assumed to be sampled on a linear radial grid starting at Rmin with spacing dR. units: Jy/sr

• Rmin (float) – Inner edge of the radial grid, i.e. the radius where the brightness is intensity[0]. units: rad

• dR (float) – Size of the cell of the radial grid, assumed linear. units: rad

• nxy (int) – Side of the square model image, which is internally computed. units: pixel

• dxy (float) – Size of the image cell, assumed equal in both x and y direction. units: rad

• u (array_like, float) – u coordinate of the visibility points where the FT has to be sampled. units: wavelength

• v (array_like, float) – v coordinate of the visibility points where the FT has to be sampled. The length of v must be equal to the length of u. units: wavelength

• vis_obs_re (array_like, float) – Real part of the observed visibilities. units: Jy

• vis_obs_im (array_like, float) – Imaginary part of the observed visibilities. The length of vis_obs_im must be equal to the length of vis_obs_re. units: Jy

• vis_obs_w (array_like, float) – Weight associated to the observed visibilities. The length of vis_obs_w must be equal to the length of vis_obs_re. units:

• dRA (float, optional) – R.A. offset w.r.t. the phase center by which the image is translated. If dRA > 0 translate the image towards the left (East). Default is 0. units: rad

• dDec (float, optional) – Dec. offset w.r.t. the phase center by which the image is translated. If dDec > 0 translate the image towards the top (North). Default is 0. units: rad

• PA (float, optional) – Position Angle, defined East of North. Default is 0. units: rad

• inc (float, optional) – Inclination of the image plane along a North-South (top-bottom) axis. If inc=0. the image is face-on; if inc=pi/2 the image is edge-on. units: rad

• check (bool, optional) – If True, check whether image and dxy satisfy Nyquist criterion for computing the synthetic visibilities in the (u, v) locations provided. Additionally check that the (u, v) points fall in the image to avoid segmentation violations. Default is False since the check might take time. For executions where speed is important, set to False.

Returns:

chi2 – The chi square, not normalized.

Return type:

float

See also

sampleProfile(), sweep()

6.3. GPU related

galario.HAVE_CUDA

Global variable (bool). It is True if the GPU libraries (galario.double_cuda and galario.single_cuda) are available, False otherwise. On a machine without a CUDA-enabled GPU it is always False.

galario.double.ngpus()

Return how many GPUs are available on the machine.

Returns:

ngpus – Number of GPUs available on the machine.

Return type:

int

galario.double.use_gpu()

Select the GPU to be used for the computation.

Typical call signature:

use_gpu(device_id)

Parameters:

device_id (int) – ID of the GPU to be used for the computation.

Notes

If more than one GPU is present, device_id might not coincide with the ID reported by the nvidia-smi command, which reflects the PCI order. If the system administrator does not provide instructions on how to set the device_id, we recommend to start from device_id=0 and simultaneously check which GPU is used with watch -n0.1 nvidia-smi.

6.4. Other useful functions

galario.double.threads()

Set and get the number of threads to be used in parallel sections of the code.

To set, pass num>0. To get the current setting, call without any argument.

Typical call signatures:

default_nthreads = threads() # in the first call, a default is preset
threads(num=16) # now change the default

Parameters:

num (int, optional) –

On the GPU, num is the square root of the number of threads per block to be used. 1D kernels are launched with linear blocks of size num*num. 2D Kernels are launched with square blocks of size num*num.

On the CPU, this sets the number of openMP threads. The default is omp_get_max_threads() which can be set through the OMP_NUM_THREADS environment variable. If compiled without openMP support, num is ignored and this function always returns 1.

Notes

The CUDA documentation suggests starting with num*num>=64 and multiples of 32, e.g. 128, 256. GPU cards with compute capability between 2 and 6.2 have maximum number of threads per block of 1024, which is achieved for num=32.

Check the maximum number of threads per block of your GPU by running the deviceQuery command.

On the CPU, it may useful to experiment with more threads than available cores to see if hyperthreading provides any benefit.

galario.double.get_image_size()

Compute the recommended image size given the (u, v) locations.

Typical call signature:

nxy, dxy = get_image_size(u, v, PB=0, f_min=5., f_max=2.5, verbose=False)

Parameters:

• u (array_like, float) – u coordinate of the visibility points where the FT has to be sampled. units: wavelength

• v (array_like, float) – v coordinate of the visibility points where the FT has to be sampled. The length of v must be equal to the length of u. units: wavelength

• PB (float, optional) – Primary beam of the antenna, e.g. 1.22*wavelength/Diameter for an idealized antenna with uniform illumination. units: rad

• f_min (float) – Size of the field of view covered by the (u, v) plane grid w.r.t. the field of view covered by the image. Recommended to be larger than 3 for better results. units: pure number

• f_max (float) – Nyquist rate: numerical factor that ensures the Nyquist criterion is satisfied when sampling the synthetic visibilities at the specified (u, v) locations. Must be larger than 2. The maximum (u, v)-distance covered is f_max times the maximum (u, v)-distance of the observed visibilities. units: pure number

• verbose (bool, optional) – If True, prints information on the criteria to be fulfilled by nxy and dxy.

Returns:

• nxy (int) – Size of the image along x and y direction. units: pixel

• dxy (float) – Returned only if not provided in input. Size of the image cell, assumed equal and uniform in both x and y direction. units: cm

galario.double.check_image_size()

Check whether the setup of the (u, v) plane satisfies Nyquist criteria for (u, v) plane sampling.

Typical call signature:

check_image_size(u, v, nxy, dxy, duv, PB=0, verbose=False)

Parameters:

• u (array_like, float) – u coordinate of the visibility points where the FT has to be sampled. units: wavelength

• v (array_like, float) – v coordinate of the visibility points where the FT has to be sampled. The length of v must be equal to the length of u. units: wavelength

• nxy (int) – Size of the image along x and y direction. units: pixel

• dxy (float) – Size of the image cell, assumed equal in both x and y direction. units: rad

• duv (float) – Size of the cell in the (u, v) plane, assumed uniform and equal on both u and v directions. units: wavelength

• verbose (bool, optional) – If True, prints information on the criteria to be fulfilled by nxy and dxy.

galario.double.get_coords_meshgrid()

Compute the (R.A, Dec.) coordinate mesh grid to create the image. (x, y) axes are the (R.A, Dec.) axes: x increases leftwards, y increases upwards. All coordinates are computed in linear pixels units. To convert to angular units, just multiply the output by the angular pixel size. To put (0, 0) coordinates offset by (Dx, Dy) w.r.t. the image center, specify Dx, Dy.

Typical call signature:

x, y, x_m, y_m, R_m = get_coords_meshgrid(nrow, ncol, dxy=dxy, inc=inc, Dx=Dx, Dy=Dy, origin='lower')

Parameters:

• nrow (int) – Number of rows of the image.

• ncol (int) – Number of columns of the image.

• dxy (float, optional) – Size of the image cell, assumed equal in both x and y direction. By default is 1, thus implying the output arrays are expressed in number of pixels. units: rad

• inc (float, optional) – Inclination along the North-South axis, default is zero. units: rad

• Dx (float, optional) – Offset of the source along the x-axis (R.A.), default is zero. If positive, moves the origin to the East (left). units: rad

• Dy (float, optional) – Offset of the source along the y-axis (Dec.), default is zero. If positive, moves the origin to the North (top). units: rad

• origin (['upper' | 'lower'], optional) – Set the [0,0] pixel index of the matrix in the upper left or lower left corner of the axes. It follows the same convention as in matplotlib matshow and imshow commands. Declination axis and the matrix y axis are parallel for origin='lower', anti-parallel for origin='upper'. The central pixel corresponding to the (RA, Dec) = (0, 0) is always [Nxy/2, Nxy/2]. For more details see the Technical Requirements page in the online docs.

Returns:

• x, y (array_like, float) – Pixel coordinates along the (R.A., Dec.) directions. units: same as dxy. If dxy=1: number of pixels.

• x_m, y_m (array_like, float) – Pixel coordinate meshgrid along the (R.A., Dec.) directions. units: same as dxy. If dxy=1: number of pixels.

• R_m (array_like, float) – Radial coordinate meshgrid. units: same as dxy. If dxy=1: number of pixels.

galario.double.sweep()

Create a 2D model image from an axisymmetric brightness profile.

The image is created assuming that the x-axis (R.A.) increases from right (West) to left (East) and the y-axis (Dec.) increases from bottom (South) to top (North).

Typical call signature:

image = sweep(intensity, Rmin, dR, nxy, dxy, inc=0)

Parameters:

• intensity (2D array_like, float) – Array containing the radial brightness profile of the model. The brightness profile is assumed to be sampled on a linear radial grid starting at Rmin and with spacing dR. units: Jy/sr

• Rmin (float) – Inner edge of the radial grid, i.e. the radius where the brightness is intensity[0]. units: rad

• dR (float) – Size of the cell of the radial grid, assumed linear. units: rad

• nxy (int) – Side of the square model image. units: pixel

• dxy (float) – Size of the image cell, assumed equal in both x and y direction. units: rad

• inc (float, optional) – Inclination of the image plane along a North-South (top-bottom) axis. If inc=0. the image is face-on; if inc=pi/2 the image is edge-on. units: rad

Returns:

image – Image of the surface brightness. units: Jy/pixel

Return type:

(nxy, nxy) array_like, float

galario.double.apply_phase_vis()

Apply phase to sampled visibility points as to translate the image in the real space by an offset dRA along Right Ascension (R.A.) and dDec along Declination. R.A. increases towards left (East), thus dRA>0 translates the image towards East.

Typical call signature:

vis_shifted = apply_phase_vis(dRA, dDec, u, v, vis)

Parameters:

• dRA (float) – Right Ascension offset. units: rad

• dDec (float) – Declination offset. units: rad

• u (array_like, float) – u-coordinates of visibility points. units: observing wavelength

• v (array_like, float) – v-coordinates of visibility points. units: observing wavelength

• vis (array_like, complex) – complex visibilities, of form Real(Vis) + i*Imag(Vis). units: Jy

Returns:

vis_out – shifted complex visibilities units: arbitrary, same as vis

Return type:

array_like, complex

galario.double.reduce_chi2()

Compute the chi square of observed and model visibilities.

Typical call signature:

chi2 = reduce_chi2(vis_obs_re, vis_obs_im, vis_obs_w, vis)

Parameters:

• vis_obs_re (array_like, float) – Real part of the observed visibilities. units: Jy

• vis_obs_im (array_like, float) – Imaginary part of the observed visibilities. The length of vis_obs_im must be equal to the length of vis_obs_re. units: Jy

• vis_obs_w (array_like, float) – Weight associated to the observed visibilities. The length of vis_obs_w must be equal to the length of vis_obs_re.

• vis (array_like, complex) – Complex model visibilities. The length of vis must be equal to the length of vis_obs_re. units: Jy

Returns:

chi2 – The chi square, not normalized.

Return type:

float

6.5. Exceptions

galario’s C++ core throws various exception. They can be distinguished by the type and the error message.

Event	C++	Python
Out of memory on GPU	std::bad_alloc	MemoryError
Invalid argument (CPU/GPU)	std::invalid_argument	ValueError
Miscellaneous, including out of memory (CPU/GPU)	std::runtime_error	RuntimeError