python - 使用HEALPy节省像素数

我有一个要绘制到healpix贴图上的星系列表（我正在使用healpy来做），每个星系都有一个固定的通量，我需要以这样的方式绘制它们，以使每个星系的通量都保存在地图上。

这是我的代码：

import numpy as np
import matplotlib.pyplot as plt
import healpy as hp
pi = np.pi

nside = 8
xsize = 100

ra = np.array([pi/4,pi/3])
dec = np.array([pi/4,pi/3])
flux = np.array([10,20])

hpm = np.zeros(hp.nside2npix(nside)) #Blank healpix map
pixindex = hp.ang2pix(nside, dec, ra)

np.add.at(hpm,pixindex,flux) #Add flux onto correct pixels

img=hp.mollview(hpm,coord=['E'],xsize=xsize,return_projected_map=True)
print(np.sum(img[img>0]))

我得到的结果是140，而不是30，这是通量的真实总和。

我知道发生了什么，并且相同的光通量分布在多个像素上（第一个星系为6个像素，第二个星系为4个像素），我知道我可以这样做：

newimg = img * (np.sum(flux)/np.sum(img[img>0]))

这样可以节省总的光子数量，但不一定节省我需要的每个星系的光子数量。也就是说，此方法最终导致第一个星系贡献了12.86的通量，第二个星系贡献了17.14的通量。

有没有一种方法可以计算出每个坐标要占用多少像素，然后根据此值改变转储的通量？

提前致谢！

最佳答案

函数xsize的hp.mollview参数应仅用于绘图目的。如果要操纵地图的分辨率，请使用hp.pixelfunc.ud_grade

例如，如果您想从nside=8转到nside=32，

下线后

np.add.at(hpm, pixindex, flux) #Add flux onto correct pixels

将ud_grade函数与power=-2一起使用，可以使总通量保持不变：

hpm_nside_32 = hp.pixelfunc.ud_grade(hpm, power=-2, nside_out=32)

总和

np.sum(hpm_nside_32)

将在30保留。

如果您需要mollview来保护通量，我无法提供解决方案。我能得到的最接近的结果是，根据mollview图像中的像素数与hpm中的像素数之比，将img值缩小。第一项xsize * xsize / 2.是Moll视图中的像素数，第二项2. * np.pi / 8.是半短轴为半长轴pi * r * 2r长度一半的椭圆的面积与面积之比矩形4r * 2r的大小。

len(hpm) / ((xsize * xsize / 2.) * (2. * np.pi / 8.)) * np.sum(img[img > 0])

当xsize = 100时，总和将成为27.380；当xsize = 1000时，近似值更好，在29.981处；当xsize = 10000时，总通量变为29.988。

一种近似的替代方法是计算-inf中非img像素数与地图中像素数（对于768为nside=8）的比率：

float(len(hpm)) / float(np.sum(img>=0)) * np.sum(img[img>0])

在xsize = 100, 1000, 10000处，通量将分别在28.295, 30.075, 29.998处。

关于python - 使用HEALPy节省像素数，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/45329503/