2016-09-15 4 views
0

J'ai mes données dans ce format datapoint[37][19] dans l'espace phi-theta. Mais parce que mes données ne peuvent pas couvrir tout le ciel, il y a du NaN dans le tableau datapoint. Il y a environ la moitié de NaN dans l'ensemble datapoint. Environ 9/10 valeurs non NaN dans datapoint sont négatifs, environ 1/10 d'entre eux sont positifs.Comment Healpix traite-t-il NaN pour interpoler des données sur la carte du ciel?

J'ai essayé cette fonction d'interpolation:

scipy.interpolate.RectSphereBivariateSpline(theta,phi,datapoint.T) 

Mais il est revenu erros. Je demande comment j'interpole les données qui contiennent NaN, les valeurs positives et négatives au niveau Healpix peut l'utiliser pour faire des cartes. J'ai une carte non lissée faite à partir de Basemap.

enter image description here

Répondre

0

Je ne sais pas si votre question est sur le point 1) traitement des données manquantes, 2) traitant de Nans ou 3) de rotation des données arbitraires sur la sphère dans une carte Healpix. 1) interpoler des données manquantes sur une grande partie du ciel nécessite au moins une certaine connaissance statistique de vos données afin de réaliser une réalisation contrainte de ce qui manque. Mais le remplissage de ces espaces n'est nécessaire que si vous effectuez des opérations non locales (comme des circonvolutions ou des dégradés), donc cela dépend de ce que vous envisagez de faire avec les données.

2) Le paramétrage des données manquantes sur NaN va certainement bousiller tous les schémas d'interpolation disponibles. 3) Le code python ci-dessous transforme un ensemble de données similaire au vôtre (pour autant que je sache) en 2 cartes Healpix, l'une utilisant un échantillonnage NGP (Nearest Grip Point), l'autre utilisant une interpolation BSpline. Notez que la seconde ne fonctionnera probablement pas en présence de NaNs, tandis que le premier est extrêmement robuste.

import healpy as hp 
import numpy as np 
import pylab as pl 

datapoint = np.zeros((37,19), dtype=np.float) 
datapoint[18,9] = 1.0 
datapoint[0,9] = -1.0 

nside  = 64 
npix  = hp.nside2npix(nside) 

# location of Healpix pixels center 
ip   = np.arange(npix) 
theta_rad, phi_rad = hp.pix2ang(nside, ip) 

# map0 : NGP sampling 
theta_deg = np.rad2deg(theta_rad) 
phi_deg = np.rad2deg(phi_rad) 
hp_0  = datapoint[np.rint(phi_deg/10.).astype(int), \ 
         np.rint(theta_deg/10.).astype(int)] 
hp.mollview(hp_0,title='NGP map') 

# map1: BSpline interpolation 
from scipy.interpolate import RectSphereBivariateSpline 
epsilon = 1.e-12 
th_in = np.linspace(epsilon, np.pi-epsilon, 19) 
ph_in = np.linspace(epsilon,2*np.pi-epsilon, 37) 
lut = RectSphereBivariateSpline(th_in, ph_in, datapoint.T, s=1) 
hp_1 = lut.ev(theta_rad, phi_rad) 
hp.mollview(hp_1,title='BSpline map') 

pl.show()