2

I need to generate a Healpyx map (using Healpy) from random $a_{\ell m}$, for a spin-2 function.

Schematically, this should look like that:

import healpy as hp
nside = 16  # for example

for el in range(1, L+1): #loop over ell mode
    for m in range(-el,el): #for each ell mode loop over m
        ind = hp.sphtfunc.Alm.getidx(nside, el, m)
        if m == 0:
            a_lm[ind] = np.random.randn()
        else:
            a_lm[ind] = np.random.randn() + 1j * np.random.randn()

a_tmp = hp.sphtfunc.alm2map(a_lm, nside, pol=True)

My two questions are:

1) how do I initialise a_lm ? Specifically, what would be its dimension, using 

a_lm = np.zeros(???)

2) if I understood correctly, the output a_tmp is a 1 dimensional list. How do I reshape it into a two-dimensional list (the map) for plotting?

johnhenry
  • 1,293
  • 5
  • 21
  • 43

2 Answers2

0

1) What properties do you want your alm to have? You could also just assume a certain power spectrum (C_ell) and use hp.synalm() or hp.synfast().

For the initialization, you've already implemented that m goes from -ell to +ell, so you have a one-dimensional array of length sum_0^ell [2ell+1]. Doing the math should give you the length you need.

2) For the plotting, you could just directly generate a random map and then use e.g. hp.mollview(), which takes the 1-dimensional HEALPix map.

Alternatively, you can use hp.alm2map() to convert your alm to a map.

I also suggest you check out the tutorial for the plotting.

Daniel Lenz
  • 3,334
  • 17
  • 36
0
  1. Usually we can follow the following steps to get the length of a_lm.
import healpy as hp

inside = 16
# Get the maximum multipole with the current nside
lmax = 3*nside - 1 #This can vary according to the use. In cosmology, the common value is 2*nside
alm_len = hp.Alm.getsize(lmax)

a_lm = np.empty(alm_len)
  1. I think the tutorial linked in @Daniel's answer is a good resource for plotting Healpix maps.
Jyothis Chandran
  • 11
  • 4