← All posts

Investigate rotation of masks in healpy

Investigating a question about healpy on Stackoverflow

import healpy as hp
import numpy as np
%matplotlib inline
hp.disable_warnings()
nside = 16
npix = hp.nside2npix(nside)
m = hp.ma(np.arange(npix, dtype=np.float32))
mask = np.zeros(npix, dtype=np.bool)
mask[hp.query_strip(nside, np.radians(75), np.radians(105))] = 1
mask[hp.query_disc(nside, hp.dir2vec(0,0, lonlat=True), np.radians(40))] = 1
m.mask = mask
hp.mollview(m);
output 1
gal2eq = hp.Rotator(coord=["G","E"])
m_rotated = gal2eq.rotate_map_pixel(m)
mask_rotated = gal2eq.rotate_map_pixel(m.mask)
hp.mollview(mask_rotated)
output 2
hp.mollview(m_rotated.mask)
output 3

Now in the first case healpy fills the map with UNSEEN and then interpolation is handled by HEALPix C++. I don’t know how internally HEALPix handles that. In the second case we pass a map of 0 and 1 and HEALPix does the interpolation, but we don’t trigger any special case of handling UNSEEN values.

# consider all values less than 1 masked
np.testing.assert_array_equal(m_rotated.mask, mask_rotated == 1)
# consider only values of 1 masked
np.testing.assert_array_equal(m_rotated.mask, mask_rotated > 0)
# try a value close to 1
np.testing.assert_array_equal(m_rotated.mask, mask_rotated > .9)
# try a value close to 1
np.testing.assert_array_equal(m_rotated.mask, mask_rotated > .999)
# try a value close to 1
np.testing.assert_array_equal(m_rotated.mask, mask_rotated > .9999)