Conserver le nombre de pixels avec HEALPy

Question

J'ai une liste de galaxies à tracer sur une carte de santé (que j'utilise healpy pour faire) chaque galaxie a un flux défini et j'ai besoin de les tracer de manière à ce que le flux pour chaque galaxie est conservée sur la carte.

Ceci est mon code:

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])) 

Le résultat est que je reçois 140 et non 30 qui est la véritable somme des flux.

je reçois ce qui se passe et que le même flux se propage sur plusieurs pixels (6 pour la première galaxie et 4 pixels pour le second) et je suis conscient que je pouvais faire:

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

et cela permettrait de conserver le nombre total de photons, mais il ne conserverait pas nécessairement le nombre de photons de chaque galaxie qui est ce dont j'ai besoin. Si cette méthode se termine avec la première galaxie contribuant un flux de 12.86 et la deuxième galaxie un flux de 17.14.

Existe-t-il un moyen de déterminer le nombre de pixels que prendra chaque coordonnée, puis de modifier la quantité de flux déposée en fonction de cela?

Merci d'avance!

Answer 1

Le paramètre xsize de la fonction hp.mollview doit être utilisé uniquement à des fins de représentation. Si vous souhaitez manipuler la résolution de la carte, utilisez hp.pixelfunc.ud_grade

Par exemple, si vous voulez aller nside=8-nside=32,

Après la ligne

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

utiliser la fonction ud_grade avec power=-2, de sorte que le flux total peut être conservé:

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

la somme

np.sum(hpm_nside_32) 

sera conservée à 30.

Je ne peux pas offrir une solution si vous avez besoin du mollview pour conserver le flux. Le plus proche que je peux obtenir est en redimensionnant la valeur img en fonction du rapport entre le nombre de pixels dans l'image mollview et le nombre de pixels dans le hpm. Le premier terme xsize * xsize/2. est le nombre de pixels dans le mollview, le second terme 2. * np.pi/8. est le rapport de l'aire d'une ellipse avec semi-mineur axe la moitié de la longueur de l'axe semi-majeur pi * r * 2r à l'aire d'un rectangle 4r * 2r.

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

Lorsque xsize = 100, la somme deviendra 27.380; quand xsize = 1000, l'approximation est beaucoup mieux, à 29.981; Lorsque xsize = 10000, le flux total devient 29.988.

Une autre façon de vous donner une bonne approximation consiste à calculer le rapport entre le nombre de non -inf pixels img au nombre de pixels dans votre carte (qui est 768 pour nside=8):

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

Au xsize = 100, 1000, 10000, le flux sera à 28.295, 30.075, 29.998 respectivement.