I have a data grid where the rows represent theta (0, pi) and the columns represent phi (0, 2*pi) and where f(theta,phi) is the density of dark matter at that location. I wanted to calculate the power spectrum for this and have decided to use healpy. What I can not understand is how to format my data for healpy to use. If someone could provide code (in python for obvious reasons) or point me to a tutorial, that would be great! I have tried my hand at doing it with the following code: #grid dimensions are Nrows*Ncols (subject to change) theta = np.linspace(0, np.pi, num=grid.shape[0])[:, None] phi = np.linspace(0, 2*np.pi, num=grid.shape[1]) nside = 512 print "Pixel area: %.2f square degrees" % hp.nside2pixarea(nside, degrees=True) pix = hp.ang2pix(nside, theta, phi) healpix_map = np.zeros(hp.nside2npix(nside), dtype=np.double) healpix_map[pix] = grid But, when I try to execute the code to do the power spectrum. Specifically, : cl = hp.anafast(healpix_map[pix], lmax=1024) I get this error: TypeError: bad number of pixels If anyone could point me to a good tutorial or help edit my code that would be great. More specifications: my data is in a 2d np array and I can change the numRows/numCols if I need to. Edit: I have solved this problem by first changing the args of anafast to healpix_map. I also improved the spacing by making my Nrows*Ncols=12*nside*nside. But, my power spectrum is still giving errors. If anyone has links to good documentation/tutorial on how to calculate the power spectrum (condition of theta/phi args), that would be incredibly helpful. Continue reading...