WcsGrid#

class cygrid.WcsGrid#

Bases: Cygrid

WCS version of Cygrid.

Parameters:

headerDictionary or anything that fits into WCS

The header must contain a valid WCS representation for a three-dimensional data cube (spatial-spatial-frequency).

dbg_messagesBoolean, optional (Default: False)

Do debugging output.

datacubearray [nD] of float32 or float64, optional (Default: None)

Provide pre-allocated or even pre-filled array for datacube.

Usually (if datacube=None, the default) a datacube object will be created automatically according to the given FITS header dictionary (for the spatial part) and the dimensions of the input raw data to be gridded (see also grid). Providing datacube manually might be worthwhile if some kind of repeated/ iterative gridding process is desired. However, for almost all use cases it will be sufficient to repeatedly call the grid method. Cygrid won’t clear the memory itself initially, so make sure to handle this correctly.

weightcubearray [nD] of float32 or float64, optional (Default: None)

As datacube but for the weight array.

dtypefloat32 or float64 (Default: float32)

Desired output format of data cube (if cygrid is allocating this; otherwise, the dtype of the input datacube provided to the constructor takes precedence).

Raises:

TypeError: Input datacube/weightcube must be numpy array. Input datacube must have floating point type. Weightcube dtype doesn’t match datacube dtype.
ShapeError: Datacube/weightcube shape doesn’t match fits header.

Examples

The following provides a minimal example. For more detailed information we refer to the user manual or the Jupyter tutorial notebooks:

from astropy.io import fits
import matplotlib.pyplot as plt
import cygrid

# read-in data
lon, lat, rawdata = get_data()

# define target FITS/WCS header
header = create_fits_header()

# prepare gridder
kernelsize_sigma = 0.2

kernel_type = 'gauss1d'
kernel_params = (kernelsize_sigma, )  # must be a tuple
kernel_support = 3 * kernelsize_sigma
hpx_maxres = kernelsize_sigma / 2

mygridder = cygrid.WcsGrid(header)
mygridder.set_kernel(
    kernel_type,
    kernel_params,
    kernel_support,
    hpx_maxres,
    )

# do the actual gridding
mygridder.grid(lon, lat, rawdata)

# query result and store to disk ...
data_cube = mygridder.get_datacube()
fits.writeto('example.fits', header=header, data=data_cube)

# ... or plot
fig = plt.figure()
ax = fig.add_subplot(111, projection=WCS(header).celestial)
ax.imshow(
    data_cube[0], origin='lower', interpolation='nearest'
    )
plt.show()

Methods Summary

`clear_cache`(self)	Clear all internal caches/dictionaries.
`clear_data_and_weights`(self)	Reset data and weights arrays.
`get_datacube`(self)	Return final data.
`get_header`(self)	Return header object for reference.
`get_pixel_coords`(self)	Return pixel coordinates of the datacube's xy-plane
`get_unweighted_datacube`(self)	Return final unweighted data.
`get_wcs`(self)	Return WCS object for reference.
`get_weights`(self)	Return final weights.
`get_world_coords`(self)	Return world coordinates of the datacube's xy-plane
`grid`(self, ndarray lons, ndarray lats, ...)	Grid irregularly positioned data points (spectra) into the data cube.
`set_kernel`(self, kernel_type, kernel_params, ...)	Set the gridding kernel type and parameters.
`set_num_threads`(self, int nthreads)	Change maximum number of threads to use.

Methods Documentation

clear_cache(self)#

Clear all internal caches/dictionaries.

Notes

To speed-up processing, cygrid uses internal caches (see Cygrid paper for details). This is usually a minor contribution to the total memory usage, but there could be scenarios, e.g., when one grids into very large maps having very small pixel sizes, that the internal cache grows too much. In this case, we recommend to sort the input data by latitude, grid the data in chunks, and call this function every now and then.

clear_data_and_weights(self)#: Reset data and weights arrays. (Caches stay filled!)

get_datacube(self)#

Return final data.

Returns:

dataarray [nD] of float32 or float64: The gridded spectral data.

Notes

The actual shape of the returned array depends on the input raw data signal shape and the target pixel array size.

get_header(self)#

Return header object for reference.

Returns:

headerdict or WCS-compatible: Header object provided to the constructor.

get_pixel_coords(self)#

Return pixel coordinates of the datacube’s xy-plane

Returns:

x/yarray [2D]: Pixel coordinates \((x, y)\) of the datacube.

get_unweighted_datacube(self)#

Return final unweighted data. (For debugging only.)

Returns:

unweighted_dataarray [nD] of float32 or float64: The gridded spectral unweighted data. The gridded data is the quotient of the unweighted data and the weights.

Notes

The actual shape of the returned array depends on the input raw data signal shape and the target pixel array size.

get_wcs(self)#

Return WCS object for reference.

Returns:

wcsWCS: WCS object created from input header dictionary.

get_weights(self)#

Return final weights.

Returns:

weightsarray [nD] of float32 or float64: The gridded spectral weights.

Notes

The actual shape of the returned array depends on the input raw data signal shape and the target pixel array size.

get_world_coords(self)#

Return world coordinates of the datacube’s xy-plane

Returns:

lons/latsarray [2D]: World coordinates \((l, b)\) of the datacube.

grid(self, ndarray lons, ndarray lats, ndarray data, ndarray weights=None)#

Grid irregularly positioned data points (spectra) into the data cube.

After successful gridding, you can obtain the resulting datacube with the get_datacube method. The associated weight cube is accessible with get_weights.

Parameters:

lons, latsarray [1D] of float64: Flat lists/arrays of input coordinates \((l, b)\).
data/weightsarray [nD] of float32 or float64: The spectra and their weights (optional) for each of the given coordinate pairs, \((l, b)\). First axis must match lons/lats size. The shape of the input data will determine the output datacube shape: the first axes of the datacube will match the ( raw) data (aka input signal) shape - only the first dimension is being stripped as it is associated with the number of coordinate pairs -, while the last axes match the shape of the target pixel array (e.g., the shape of the target map).

Raises:

ShapeError: Input coordinates/data points length mismatch. Number of spectral channels mismatch.

Notes

Internally, cygrid always works with a 3D representation of the data cube, i.e., with (possibly redundant) spectral axis. However, as the original shape is stored, the gridded data returned by get_datacube will always be consistent with the input raw data signal (in terms of dimensions).
All input parameters need to be C-contiguous (cygrid will re-cast if necessary).

set_kernel(self, kernel_type, kernel_params, double sphere_radius, double hpx_max_resolution)#

Set the gridding kernel type and parameters.

Parameters:

kernel_typestr

Set the kernel type.

The following names/types are available: ‘gauss1d’, ‘gauss2d’, ‘tapered_sinc’ (see Notes for details)

kernel_paramsarray

Set the kernel parameters for the chosen type (see Notes for details)

sphere_radiusdouble

Kernel sphere radius.

This is controls out to which distance the kernel is computed for. For Gaussian kernels, values much larger than \(3\ldots4 \sigma_\mathrm{k}\) do not make much sense.

hpx_max_resolutiondouble

Maximum acceptable HPX resolution (\(\sigma_\mathrm{k}^{[\mathrm{maj}]} / 2\) is a reasonable value).

Notes

Below you find a list of kernel-names and required parameters:

'gauss1d', (kernel_sigma,)
'gauss2d', (kernel_sigma_maj, kernel_sigma_min, PA)
'tapered_sinc', (kernel_sigma, param_a, param_b)

All numbers are in units of degrees, except for PA, param_a and param_b. PA (the position angle) is in units of radians (for efficiency). param_a and param_b should be 2.52 and 1.55, respectively, for optimal results!

The kernel size (sigma) defines the amount of “smoothness” applied to the data. If in doubt a good value is about 50% of the true/input angular resolution of the data (this will result in about 10% degradation of the final angular resolution.)

set_num_threads(self, int nthreads)#

Change maximum number of threads to use.

This is a convenience function, to call omp_set_num_threads(), which is otherwise not possible during runtime from Python.