llama.files.healpix.utils module

Utility functions used across healpix_skymap classes.

llama.files.healpix.utils.check_valid_nside(nside)

Checks whether nside is a valid HEALPix NSIDE value. Raises a ValueError if it is not.

Parameters

nside (int or array) – An integer or array of integers representing HEALPix NSIDE values.

Raises

ValueError – If any of the values are not valid HEALPix NSIDE values.

llama.files.healpix.utils.check_valid_nuniq(indices)

Checks that indices are valid NUNIQ indices.

Raises

ValueError – If indices are not valid NUNIQ indices, i.e. they are not integers greater than 3.

llama.files.healpix.utils.dangle_rad(ra, dec, mapra, mapdec)

Get an array of angular distances in radians between the point specified in ra, dec and the points specified in mapra, mapdec. All angles, including the return value, are specified in radians.

sin(DEC1)*sin(DEC2) + cos(DEC1)*cos(DEC2)*cos(RA2-RA1)

So we get the angle by taking the arccos of this function.

Parameters
  • ra (float) – Right ascension of the single point [radians]

  • dec (float) – Declination of the single point [radians]

  • mapra (array) – Right ascensions whose angular distances to the single point will be calculated; must be same length as mapdec. [radians]

  • mapdec (array) – Declinations whose angular distances to the single point will be calculated; must be same length as mapra. [radians]

Returns

dangles – Angular distances between the single point provided and the arrays of points [radians].

Return type

array

llama.files.healpix.utils.downres_mask_nest(mask_nested, mask_nside, nside, mask_sorted=False)

Get a lower-resolution set of pixel indices that fully cover the region specified in mask_nested. The returned indices correspond to larger pixels which might contain part of the sky that are not included in the original mask (e.g. if the larger pixel only contains one of the original pixels in mask_nested).

Parameters
  • mask_nested (array-like) – The set of nested HEALPix pixels that you’d like to have at lower resolution.

  • mask_nside (int) – The HEALPix nside parameter of this mask. Must be a power of 2.

  • nside (int) – The nside value you would like for your output. Must be a power of 2.

  • mask_sorted (bool, optional) – Whether the input mask_nested is already sorted. If it is sorted and is large, specifying this will skip a sort step, potentially saving time.

Returns

lowres_mask_nested – The smallest sky region that can be specified at the resolution given by nside that contains the sky region specified in mask_nested at mask_nside. This will be roughly a factor of (nside/mask_nside)**2 smaller, but do not count on this ratio, since the pixels in the returned mask might not have all of their corresponding subpixels in the original mask_nested; in pathological cases, the returned array can be up to the same size (if every pixel in mask_nested only appears in a single larger pixel in at resolution nside). Will be a new array independent of inputs, so feel free to mutate it.

Return type

np.array

Raises

ValueError – if mask_nside is smaller than nside

llama.files.healpix.utils.nest2uniq(indices, nside)

Return the NUNIQ pixel indices for nested indices with NSIDE=``nside``.

Parameters
  • indices (array) – Indices of HEALPix pixels in NEST ordering.

  • nside (int) – A valid NSIDE value not greater than any of the NSIDE values of the indices.

Returns

uniq – An array of HEALPix pixels in NEST ordering corresponding to the input indices and nside.

Return type

array

Raises

ValueError – If nside is an invalid NSIDE value, i.e. not a power of 2.

Examples

>>> import numpy as np
>>> nest2uniq(np.array([284, 286, 287,  10,   8,   2, 279, 278, 285]), 8)
array([540, 542, 543, 266, 264, 258, 535, 534, 541])
llama.files.healpix.utils.read_partial_skymap(infile, mask_nested, mask_nside, mask_sorted=False, memmap=True)

Read in pixels from a FITS skymap (or a gzip-compressed FITS skymap) that lie in a specific sky region. Attempts to minimize memory usage by memory-mapping pixels in the input file and only loading those specified in mask_nested.

Parameters
  • infile (str) – A FITS HEALPix skymap (optionally compressed; see note under memmap)

  • mask_nested (bool, optional) – HEALPix pixel indices in nested ordering specifying the region of the skymap that should be loaded

  • mask_nside (int) – The resolution of the pixel mask

  • mask_nested – Whether mask_nested is already sorted. You can specify mask_sorted=True, potentially saving some sorting time

  • memmap (bool, optional) – If True, use memory mapping when reading in files. This can VASTLY reduce memory required for high-resolution skymaps. If infile is gzip-compressed and memmap is True, then infile will be decompressed to a temporary file and data will be read from it (necessary to constrain memory usage); for high-resolution skymaps, this can require the availability of several gigabytes of tempfile storage.

Returns

partial_skymap – A partial skymap table in nested ordering. Has two columns: INDICES and PROB. If the resolution of the original HEALPix skymap file is lower than that of the mask, then any pixels overlapping with those in mask_nested will be used; this might result in a larger portion of the skymap being used than that specified in mask_nested. The resolution of this skymap will be the resolution of the smallest pixel loaded from the input file (in the case of ring or nested ordering, this is just the resolution of the input skymap).

Return type

astropy.table.Table

llama.files.healpix.utils.sky_area()

Get the area of the entire sky as an Astropy.units.Quantity.

llama.files.healpix.utils.sky_area_deg()

Get the area of the entire sky in square degrees.

llama.files.healpix.utils.uniq2nest_and_nside(indices)
Parameters

indices (array) – A scalar or numpy array of NUNIQ indices

Returns

  • nest_indices (array) – The indices in nest ordering at their respective resolutions

  • nside (array) – The resolution expressed as NSIDE of the indices given in nest_indices

Examples

>>> import numpy as np
>>> nuniq = np.array([540, 542, 543, 266, 264, 258, 535, 534, 541])
>>> uniq2nest_and_nside(nuniq)
(array([284, 286, 287,  10,   8,   2, 279, 278, 285]), array([8, 8, 8, 8, 8, 8, 8, 8, 8]))

Confirm that this is the inverse of nest2uniq: >>> all(nest2uniq(*uniq2nest_and_nside(nuniq)) == nuniq) True

llama.files.healpix.utils.uniq2nside(indices)

Get the NSIDE value of the given NUNIQ-ordered indices.

Raises

ValueError – If indices are not valid NUNIQ indices, i.e. they are not integers greater than 3.

llama.files.healpix.utils.uniq2uniq_lowres(indices, nside)

Find the indices at resolutions nside (each of which MUST be lower resolution than that of the corresponding pixel in indices) which contain the higher resolution pixels described by the corresponding indices.

Parameters
  • indices (array) – Indices of HEALPix pixels in NUNIQ ordering.

  • nside (array) – A valid NSIDE value not greater than any of the NSIDE values of the indices. Each of these corresponds to an index in indices.

Returns

lowres_indices – Indices of HEALPix pixels in NUNIQ ordering with resolution nside which contain the corresponding pixels in indices.

Return type

array

Raises

ValueError – If the NSIDE resolution of any index in indices is smaller than nside, or if nside is an invalid NSIDE value, i.e. not a power of 2.

Examples

Take a few NUNIQ pixels at NSIDE=16 and an NUNIQ pixel at NSIDE=1024 and find which NUNIQ pixels they correspond to at NSIDE=8. The first four should all correspond to the same pixel at NSIDE=8.

>>> import numpy as np
>>> indices = np.array([1024, 1025, 1026, 11295603])
>>> uniq2uniq_lowres(indices, 8)
array([256, 256, 256, 689])
llama.files.healpix.utils.uniq_intersection(indices1, indices2)

Downselect the pixel indices given in indices1 to the set that overlaps with pixels in indices2 and return pairs of indices into both of the input index lists showing which pixel in indices2 each downselected pixel from indices1 overlaps with. Use this to get rid of pixels outside of a given region, or alternatively use it as part of a calculation involving two multi-resolution skymaps whose pixel sizes are non-identical. Should have stable O(len(indices1)*len(indices2)) performance, so make sure that your index lists are not both huge.

Parameters
  • indices1 (array) – Indices of HEALPix pixels in NUNIQ ordering.

  • indices2 (array) – Indices of HEALPix pixels in NUNIQ ordering.

Returns

  • indices_into_indices1 (array) – Indices into indices1 that overlap with indices2.

  • indices_into_indices2 (array) – Corresponding indices into indices2 that overlap with indices1.

Examples

Some pixels with NSIDE of 16, 32, and 64, respectively: >>> import numpy as np >>> indices1 = np.array([1024, 4100, 44096])

Pixels at NSIDE = 32 that overlap with only the first and last pixels of indices1: >>> indices2 = np.array([4096, 4097, 11024])

We should see correspondence between index 0 of indices1 and indices 0, 1 of indices1; and correspondence between index 2 of indices1 and index 2 of indices2: >>> uniq_intersection(indices1, indices2) (array([0, 0, 2]), array([0, 1, 2]))

llama.files.healpix.utils.uniq_rasterize(indices, values, density=True, nside=None)

Increase resolution of an NUNIQ HEALPix skymap and return NUNIQ indices and values at the same higher resolution (same NSIDE). Doesn’t do any interpolation; just used to make calculations easier.

Parameters
  • indices (array) – Non-overlapping HEALPix NUNIQ indices into the skymap.

  • values (array) – Pixel values corresponding to the indices.

  • density (bool, optional) – Whether the pixel values are some sort of density. If True (the default), the pixel values will not be modified (since they don’t depend on pixel area. Otherwise, the pixel value will be split evenly between sub-pixels.

  • nside (int, optional) – The NSIDE value to convert the skymap to. If not provided, the skymap will be converted to the largest NSIDE of the provided pixels.

Returns

  • raster_indices (array) – A unique flat array of HEALPix NUNIQ indices into the rasterized skymap, sorted in ascending order.

  • raster_values (array) – A flat array of pixel values of the rasterized skymap corresponding to raster_indices.

Raises

ValueError – If the given nside value is either an invalid NSIDE value or if it is lower than the highest NSIDE value calculated from indices, or if the given indices are overlapping.

Examples

Some example inputs with NSIDE of 16, 32, and 64, respectively: >>> import numpy as np >>> indices = np.array([1024, 4100, 44096]) >>> values = np.ones_like(indices)

Rasterize to the max resolution of NSIDE=64, assuming values represent a density: >>> uniq_rasterize(indices, values) (array([16384, 16385, 16386, 16387, 16388, 16389, 16390, 16391, 16392,

16393, 16394, 16395, 16396, 16397, 16398, 16399, 16400, 16401, 16402, 16403, 44096]), array([1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]))

Rasterize just the first two indices, this time explicitly to NSIDE=64 (vs. the default max NSIDE, 32 in this case), explicitly treating the values NOT as a density (and therefore splitting them amongst subpixels): >>> uniq_rasterize(indices[0:2], values[0:2], density=False, nside=64) (array([16384, 16385, 16386, 16387, 16388, 16389, 16390, 16391, 16392,

16393, 16394, 16395, 16396, 16397, 16398, 16399, 16400, 16401, 16402, 16403]), array([0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.0625, 0.25 , 0.25 , 0.25 , 0.25 ]))

llama.files.healpix.utils.upres_nest(mask_nested, mask_nside, nside, mask_sorted=False)

Get nested indices for a skymap mask at a higher resolution. Output will be sorted.

Parameters
  • mask_nested (array-like) – The set of nested HEALPix pixels that you’d like to have at higher resolution.

  • mask_nside (int) – The HEALPix nside parameter of this mask. Must be a power of 2.

  • nside (int) – The nside value you would like for your output. Must be a power of 2.

  • mask_sorted (bool, optional) – Whether the input mask_nested is already sorted. If it is sorted and is large, specifying this will skip a sort step, potentially saving time.

Returns

highres_mask_nested – The same sky region specified in mask_nested at mask_nside but at the resolution specified in nside. This will be a factor of (nside/mask_nside)**2 larger. Will be a new array independent of inputs, so feel free to mutate it.

Return type

np.array

Raises

ValueError – if mask_nside is larger than nside