Построение массива numpy в healpy

Я пытаюсь создать луч на карте healpix, используя healpy. Для начала я хотел бы иметь возможность создавать 2D-гауссовскую проекцию Моллвейда, но я действительно не знаю, с чего начать.

Я могу определить двумерный гаусс:

import numpy as np
def gaussian_2D(x,y,mu_x=0.,mu_y=0.,sig_x=1.,sig_y=1.):
    return np.exp(-0.5*(((x-mu_x) / sig_x)**2 + ((y-mu_y) / sig_y)**2))

так что я могу создать трехмерное пространство X, Y, Z, например:

delta = 0.025
x = np.arange(-4, 4, delta)
y = np.arange(-4, 4, delta)
X, Y = np.meshgrid(x,y)
Z = gaussian_2D(X,Y)

но отсюда я довольно потерян и не могу отследить много полезной документации о том, как и/или что проектировать. Любые предложения по направлению атаки будут высоко оценены!


python numpy arrays astronomy healpy
person Saul Aryeh Kohn    schedule 28.08.2014    source источник
comment
healpy использует пикселизацию HEALPix, поэтому карта представляет собой одномерный массив, где индексы соответствуют пикселям. Если вам просто нужна проекция Моллвейде, вы можете просто использовать matplotlib, см. matplotlib.org/examples/pylab_examples/ geo_demo.html   -  person Andrea Zonca    schedule 29.08.2014

Ответы (1)


arrow_upward
0
arrow_downward

вот как я это делаю:

используя небольшую хитрость. Я вставляю точку в желаемый гауссовский центр, а затем использую «размазывание», чтобы создать гауссиан с некоторой сигмой.

Вот пример:

#!/usr/bin/env python
import numpy as np
import healpy as hp
import pylab as pl

NSIDE=512 #the map garannularity

m_sm=np.arange(hp.nside2npix(NSIDE)) # creates the map
m_sm=m_sm*0. # sets all values to zero

theta=np.radians(80.) # coordinates for the gaussian
phi=np.radians(20.)

indx=hp.pixelfunc.ang2pix(NSIDE,theta,phi) # getting the index of the point corresponding to the coordinates
m_sm[indx]=1. # setting that point value to 1.

gmap=hp.smoothing(m_sm, sigma=np.radians(20.),verbose=False,lmax=1024) # creating a new map, smmeared version of m_sm

hp.mollview(gmap, title="Gaussian Map") #draw it
pl.show()

теперь, если вы хотите сделать это вручную, вы должны использовать функцию для гауссова

1) вы скармливаете ему какие-то координаты

2) вы получаете индекс, соответствующий этой координате, используя:

indx=hp.pixelfunc.ang2pix(NSIDE,theta,phi)

3) вы устанавливаете значение для этой точки на значение из вашей функции Гаусса. то есть:

my_healpy_map[indx]=my_gauss(theta, phy, mean_theta, mean_phy, sigma_theta, sigma_phy)
person Asen Christov    schedule 11.12.2014