pythonsmoothinghealpy

How to efficiently perform a top-hat (disk-like) smoothing on a healpix map?

I've got a high-resolution healpix map (nside = 4096) that I want to smooth in disks of a given radius, let's say 10 arcmin.

Being very new to healpy and having read the documentation I found that one - not so good - way to do this was to perform a "cone search", that is to find around each pixels the ones inside the disk, average them and give this new value to the pixel at the center. However this is very time-consuming.

import numpy as np
import healpy as hp

kappa = hp.read_map("zs_1.0334.fits") #Reading my file

NSIDE = 4096

t = 0.00290888  #10 arcmin
new_array = []
n = len(kappa)
for i in range(n):
     a = hp.query_disc(NSIDE,hp.pix2vec(NSIDE,i),t)
     new_array.append(np.mean(kappa[a]))

I think the healpy.sphtfunc.smoothing function could be of some help as it states that you can enter any custom beam window function but I don't understand how this works at all...

Thanks a lot for your help !

Solution

  • As suggested, I can easily make use of the healpy.sphtfunc.smoothing function by specifying a custom (circular) beam window.

    To compute the beam window, which was my problem, healpy.sphtfunc.beam2bl is very useful and simple in the case of a top-hat.

    The appropriated l_max would roughly be 2*Nside but it can be smaller depending on specific maps. One could for example compute the angular power-spectra (the Cls) and check if it dampens for smaller l than l_max which could help gain some more time.

    Thanks a lot to everyone who helped in the comments section!