python - 使用healpy进行角度分箱

Question

假设我有三个量：theta、phi 和 v(theta,phi)。我想使用角度分箱，以便我可以插入任何未来的 theta 和 phi 来获得 v。我对 healpix 完全陌生，不明白如何去做。本质上，我想要一个 theta 和 phi 的网格，然后想使用 scipy.griddata 进行插值。谢谢。

Answer 1

您可以只使用插值scipy.interpolate.interp2d，请参阅https://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.interpolate.interp2d.html，甚至不使用healpy.

但是，让我展示你如何使用healpy地图来做到这一点。这个想法是，您预先计算v(theta,phi)地图的每个像素，然后为将来计算它们属于哪个像素，theta并phi使用 非常快速地获取该像素的地图值healpy。

在这里查看我的笔记本：https ://gist.github.com/zonca/e16fcf42e23e30fb2bc7301482f4683f

我复制了下面的代码以供参考：

import healpy as hp
import numpy as np

NSIDE = 64

print("Angular resolution is {:.2f} arcmin".format(hp.nside2resol(NSIDE, arcmin=True)))
NPIX = hp.nside2npix(NSIDE)
print("Number of pixels", NPIX)
pixel_indices = np.arange(NPIX)
theta_pix, phi_pix = hp.pix2ang(NSIDE, pixel_indices)
def v(theta, phi):
    return theta ** 2
v_map = v(theta_pix, phi_pix)
new_theta, new_phi = np.radians(30), np.radians(0)
new_pix = hp.ang2pix(NSIDE, new_theta, new_phi)
print("New theta and phi are hitting pixel", new_pix)
# We find what pixel the new coordinates are hitting and get the precomputed value of V at that pixel, to increase accuracy, you can increase NSIDE of the precomputed map.
v_map[new_pix]
v(new_theta, new_phi)